arXiv:1505.06209v2 [astro-ph.GA] 28 May 2015 


Draft version May 29, 2015 

Preprint typeset using DTeX style emulateapj v. Qbll2llA 


ULTRAFAINT DWARF GALAXIES - THE LOWEST MASS RELICS EROM BEEORE REIONIZATION 


Joss Bland-Hawthorn 

Sydney Institute for Astronomy, School of Physics, University of Sydney, NSW 2006, Australia 


Ralph Sutherland 

Research School of Astronomy &: Astrophysics, Australian National University, Cotter Rd, Weston, ACT 2611, Australia 


David Webster 

Sydney Institute for Astronomy, School of Physics, University of Sydney, NSW 2006, Australia 

Draft version May 29, 2015 

ABSTRACT 

New observations suggest that ultrafaint dwarf galaxies (UED) - the least luminous systems bound 
by dark matter halos (< 10^ Lq) - may have formed before reionization. The extrapolated virial 
masses today are uncertain with estimates ranging from 10^ Mq to as high as 10^ Mq depending on 
the assumed form of the underlying potential. We show that the progenitor halo masses of UEDs can, 
in principle, be as low as Mvir ~ 10^ Mq. Under the right conditions, such a halo can survive the 
energy input of a supernova (SN) and its radiative progenitor. A clumpy (fractal) medium is much 
less susceptible to both internal and external injections of energy. It is less prone to SN sweeping 
(particularly if it is off-centred) because the coupling efficiency of the explosive energy is much lower 
than for a diffuse interstellar medium. With the aid of the 3D hydro/ionization code Fyris, we 
show that sufficient baryons are retained to form stars following a single supernova event in dark 
matter halos down to Mvir ~ 10^ Mq in the presence of radiative cooling. In these models, the gas 
survives the SN explosion, is enriched with the specific abundance yields of the discrete events, and 
reaches surface densities where low mass stars can form. Our highest resolution simulations reveal 
why cooling is so effective in retaining gas compared to any other factor. In the early stages, the 
super-hot metal-enriched SN ejecta exhibit strong cooling, leading to much of the explosive energy 
being lost. Consistent with earlier work, the baryons do not survive in smooth or adiabatic models in 
the event of a supernova. The smallest galaxies may not contribute a large fraction of matter to the 
formation of galaxies, but they carry signatures of the earliest epochs of star formation, as we show. 
These signatures may allow us to distinguish a small primordial galaxy from one that was stripped 
down to its present size through tidal interaction. We discuss these results in the context of local 
ultra-faint dwarfs and damped Lya systems {z ^ 2) at very low metallicity ([Ee/H] ^ —3). We show 
that both classes of objects are consistent with primordial low-mass systems that have experienced 
only a few enrichment events. 

Subject headings: dwarf galaxies - stellar populations - star clusters - elemental abundances 


I. INTRODUCTION 

One of the most important questions in astrophysics 
is what constituted the earliest baryonic systems, and 
whether relics of these survive to the present day. Pri¬ 
mordial objects are important to identify at any redshift 
because these retain chemical signatures of the first and 
second generations of stars (Karlsson et al. 2013, here¬ 
after KI3). Eurthermore, their existence and distribu¬ 
tion provide important constraints on warm dark mat¬ 
ter models (Kennedy et al. 2014) and how reionization 
propagated through the local universe (Busha et al. 2010; 
Ocvirk et al. 2013). 

The earliest baryonic systems must have contained 
dark matter to allow gas cooling to proceed in order to 
form the first star (Abel et al. 2002). Gas cooling and 
accretion may have started as early as z ^ 30 in a dark 
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matter halo with a mass of roughly 10^“^ Mq. The first 
stars are thought to have been very massive, short lived, 
and unique to their time (Yoshida et al. 2006; Gao et al. 
2007). There is an emerging consensus on the sequence 
of events involved in forming the first stars because the 
hydrodynamic processes appear to be relatively well de¬ 
fined (Bromm & Yoshida 2011). The same cannot be 
said about their chemical yields (K13). Some of these 
stars may have seeded the first black holes in their cores. 
An unknown combination of metals is expected to fall 
back onto the black hole seed (Podsiadlowski et al. 2002; 
Umeda & Nomoto 2003), thus complicating any attempt 
to deduce the first stellar yields. 

While many questions remain concerning the forma¬ 
tion of the first stars (Clark et al. 2011; Greif et al. 2011; 
O’Shea & Norman 2007; Norman 2010), the processes 
leading to ‘second stars’ and subsequent generations are 
even less well understood (Ritter et al. 2012; Jeon et al. 
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2014). Whether the first stars inhibit or aid the forma¬ 
tion of the second generation is a complex problem (Ahn 
& Shapiro 2007). In particular, the emergence of long- 
lived, low mass stars depends in part on whether dust 
formed from the ejected material of the first stars (Clark 
et al. 2008). If dust did not form in appreciable amounts, 
the cosmic microwave background inhibits gas cooling, 
in which case low mass stars may not have emerged un¬ 
til z < 10 when mean metallicities reached ^ 10“^ Zq 
(Schneider & Omukai 2010; Bromm & Yoshida 2011). 

The first baryonic systems. So what were the first bary- 
onic systems?^ Were these exclusively the domain of the 
most massive stars, or did they include low to interme¬ 
diate mass stars (Tsuribe & Omukai 2008; Clark et al. 
2008), and do any of these stars survive to the present 
day (Okrochkov & Tumlinson 2010)? There are few if 
any reliable observational constraints thus far on the ini¬ 
tial mass function (IMF) of the earliest galaxies. In this 
paper, our goal is to explore astrophysical situations in 
which very low mass galaxies, i.e. baryons confined by 
dark matter, have survived from ancient times to the 
present day. These objects are of great interest because 
they must have preserved chemical signatures of early 
star formation and therefore could provide direct con¬ 
straints on the properties of the earliest generations of 
stars (K13). 

One of the most interesting developments in recent 
years is the discovery of ultra-faint dwarf galaxies in 
the Galactic halo (Willman et al. 2005; Belokurov et al. 
2006; Zucker et al. 2006; Belokurov et al. 2007; Walsh 
et al. 2007). These dark-matter dominated objects have 
very low baryonic mass (< 10^ Mq), are exceedingly faint 
(< 10^ Lq), and have stellar metallicities below [Fe/H]=- 
2 (Martin et al. 2008). Remarkably, some UFDs are 
less luminous than globular clusters, but the former have 
dark matter whereas the latter do not, and hence must 
be considered as galaxies. The characteristic halo mass 
of UFDs today is highly uncertain with estimates rang¬ 
ing from 10^ to 10^ Mq (Simon & Geha 2007; Strigari 
et al. 2008; Martin et al. 2008; Bovill & Ricotti 2009; 
Wolf et al. 2010; Frebel et al. 2014). 

The reionization epoch. UV radiation from a dispersed 
population of ‘first stars’ can dissociate H 2 , the primary 
coolant responsible for cooling within dark haloes. The 
diffuse ionised gas is then too warm to stay bound to low- 
mass haloes. In order for gas to be confined, the dark- 
halo mass must have exceeded Merit ~ lO^’^ Mq if the 
gas is to cool and eventually form stars after reionization 
has taken hold (Rees 1986; Ikeuchi 1986; Efstathiou 1992; 
Barkana & Loeb 1999). These simple models based on 
ionization balance and heating are supported by 3D nu¬ 
merical simulations (Okamoto et al. 2008). Thus know¬ 
ing the characteristic halo mass of UFDs - at the time of 
formation - is important for ascertaining whether they 
formed before or after the onset of reionization. 

A more nuanced approach to reionization considers 
not only the impact of positive feedback, but also the 

^ How globular clusters fit into the early history of galaxy forma¬ 
tion, and what they may tell us about the first stellar generations, 
remains a mystery. These ancient systems, which are relatively 
metal-enriched ([Fe/H]>-2.5) and free of dark matter, were long 
considered to be the oldest baryonic relics of the early universe 
(Peebles &: Dicke 1968). 


local environment of the dwarf population. In particu¬ 
lar, reionization must have been very patchy in the lo¬ 
cal universe due to the clumpy distribution of the Local 
Group, nearby groups and the Virgo cluster. Just how 
the expanding Stromgren spheres from the early galaxies 
evolved and overlapped is difficult to establish and pre¬ 
sumably awaits future observations of the distribution 
of neutral gas at high redshift. Numerical simulations 
reveal that local reionization may have started around 
z = 11 — 14 and have completed hy z = 6 — 8 (Busha 
et al. 2010; Oevirk et al. 2013), a delay of order 400 Myr. 
A major uncertainty is the impact of the ionization front 
from Virgo: if the first stars formed significantly earlier, 
the local reionization epoch may have been earlier still. 
Regardless, much of the dwarf population (and presum¬ 
ably the globular cluster population) may have acquired 
most of their stars before reionization was complete (Ri¬ 
cotti & Gnedin 2005; Busha et al. 2010). Indeed, their 
star formation may have been curtailed by the rising lo¬ 
cal UV background at late times. 

The origin of UFDs. The nature of UFDs, whether they 
precede or follow the epoch of reionization, is an area of 
active research (Brown et al. 2014; Weisz et al. 2014). It 
is difficult to relate low-mass galaxies today with their 
early universe counterparts. Some UFDs may be the 
dynamically stripped cores of more massive systems al¬ 
though the majority, as presented here, may be relics that 
have preserved their halo mass over billions of years, with 
significant accretion of dark matter at later times. (It is 
conceivable that the inner halo UFDs have lost the mat¬ 
ter they accreted after z = 10 through tidal stripping, 
such that the halo masses inferred today broadly reflect 
their masses at the time of formation.) 

The mass of the proto-Galaxy at z = 10 would have 
been of roughly 5 x 10^ Mq with a virial radius of about 
5 kpc. This early system was surrounded by hundreds 
of subhalos with masses < 10^ Mq. Over cosmic time, a 
small fraction of halos on highly radial orbits would have 
been disrupted and absorbed by the evolving Galaxy. 
Their orbit distribution and mass ratio compared to the 
progenitor halo imply that dynamical friction is rela¬ 
tively unimportant, thus many are likely to survive to the 
present day. Some fraction of these would have accreted 
gas (Susa & Umemura 2004a,b) although the details 
will depend on many environmental factors (Benftez- 
Llambay et al. 2013). It is likely that the present popu¬ 
lation of a few dozen UFDs is a small fraction of a much 
bigger population to be revealed in future deep photo¬ 
metric surveys (Bullock et al. 2010; Hargis et al. 2014). 

In an important new study of six UFDs (Brown et al. 
2014), five objects are found to have formed 75% of their 
stellar mass before z 10 (13.2 Gyr ago) and all of their 
stars by z ^ 3 (II.6 Gyr ago). There is a strong sugges¬ 
tion here of rapid quenching by a photoionising event: it 
is plausible that these UFDs preceded the main epoch 
of local or possibly global reionization (see also Webster 
et al. 2015a). In the work that follows, this is our starting 
assumption. We perform hydro dynamical simulations of 
early progenitors of UFDs with the aim of deriving their 
dynamical, photometric and chemical properties today. 
We also consider whether our models are able to explain 
the new class of very metal poor, damped Lya systems 
at z 2 which have been argued to be intrinsically low- 
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mass systems (Cooke et al. 2011b, 2015), although rather 
less is known about these intriguing systems. Where our 
work differs from earlier research is that we focus on the 
retention of baryons — rather than disruption — in dark 
matter halos that are traditionally considered to be too 
small to retain baryons. This requires most of the action 
to take place before reionization, and in clumpy media 
that are less susceptible to supernova feedback. 

New models. In our new low-mass halo models presented 
below, we consider the recent claim that star formation in 
UFDs finished before reionization (Brown et al. 2014) al¬ 
though the isochrone ages are sufficiently uncertain that 
some of it may have continued long after (Weisz et al. 
2014). For star formation to have started so soon, some 
gas must have cooled and settled into the potential well 
before the onset of reionization. We stress that complete 
reionization does not eompletely exelude ongoing star for¬ 
mation for masses below Merit (e.g. Susa & Umemura 
2004a,b; Dijkstra et al. 2004). For it to continue to later 
times, the nuclear gas density must be high enough to 
shield the central regions from a global ionizing event. 

Consider the time it takes for the gas to evaporate from 
the low mass halos. For a cosmic ionizing UV intensity 
Jo = 10“^^ erg cm“^ s“^ Hz“^ sr“^, with a moderate 
power-law form of fy oc we obtain an ionization fiux 
of ^ 2 X 10^ photons s“^ cm“^. The rate of evaporation 
of hydrogen atoms off the surface of a confined gas cloud, 
of radius 50 pc (^ 3 x 10^^ cm^), assuming one ionising 
photon evaporates one hydrogen atom, is then given by 

mevap ~ 1-5 X 10“^(nr + 1)“^ Mq yr“\ (1) 

where the term allows for the number of recombi¬ 
nations in the escaping wind. Thus, the timescale for 
evaporation for the Mb ^ 10^ Mq of gas in the models 
shown here (revap ^ AMb/rhevap) can exceed 100 Myr, 
which is long enough for star formation to proceed in the 
interior while the outer regions are evaporating. With 
reference to positive feedback processes, various authors 
have noted that H 2 formation is enhaneed in the limit 
of a weak ionizing intensity reaching the dense gas core 
(Kang & Shapiro 1992; Tajiri & Umemura 1998; Susa & 
Umemura 2000, 2004a) which can further accelerate star 
formation. Some of these processes may assist ongoing 
star formation in UFDs after the onset of reionization 
(Brown et al. 2014; Weisz et al. 2014). 

For dwarf galaxies in the outer halo at the time of 
their formation, their internal radiation fields are likely 
to have dominated local reionization (Oevirk & Aubert 
2011). Here we consider the impact of UV radiation from 
massive young stars born within the low mass halo, and 
the impact of the subsequent supernova explosion on gas 
retention. An alternative route to suppressing star for¬ 
mation in low mass systems is to remove or heat up the 
gas through starburst-driven winds (Mac Low & Ferrara 
1999). In their 2D axisymmetric models, constant energy 
input corresponding to central explosions is found to cou¬ 
ple efficiently to a smooth gas component confined by a 
fiattened disk. Once again, these models rule out any 
gas surviving in a dark-matter halo with a mass < Merit • 
But these authors did not consider discrete SN events and 
modelled the winds moving through a smooth medium 
rather than a more physical clumpy medium where the 
hot wind fiuid couples less efficiently with the dense in¬ 


terstellar medium (Cooper et al. 2008). We consider the 
impact of winds in smooth and clumpy media, and winds 
that are off-centred within the galaxy. Both off-centred 
explosions and clumpy media are likely at any epoch and 
therefore both cases must be considered in detail. 

Here we investigate the conditions under which a ‘min¬ 
imum mass’ system confined by a dark-matter halo, 
specifically a relic of the early universe, survives to the 
present day. We consider feedback arising from high- 
mass stars due to their stellar winds, their radiation 
fields, and their subsequent supernova explosions. We 
find that these limiting cosmological systems may have 
had total masses of Mmin ^ 10^ Mq, well below the 
limit of what is thought to have survived the reionization 
epoch. In our follow-up papers, we perform new calcula¬ 
tions of chemical evolution within these ‘minimum mass’ 
galaxies over cosmic time (Webster et al. 2014, 2015a,b). 
Conceivably, one can hope to learn about the star for¬ 
mation before, during and immediately after the reion¬ 
ization epoch from distinct chemical signatures in these 
low mass systems (e.g. Bland-Hawthorn et al. 2010). 

In our new models, the initial baryon fraction is 
roughly ^ 10% of the virial mass although the final 
baryon fraction can be orders of magnitude smaller. This 
is consistent with observations of dark-matter dominated 
systems observed today in the low mass limit (McGaugh 
et al. 2010). As we shall see, our new models explain 
some ultra-faint dwarf galaxies (Willman et al. 2005; Ir¬ 
win et al. 2007; Simon & Geha 2007; Martin et al. 2008; 
Strigari et al. 2008) and the peculiar chemical signatures 
observed in very metal poor, damped Lya systems (Pet- 
tini et al. 2008; Penprase et al. 2010; Cooke et al. 2011a). 
We discuss our new simulations in light of both classes 
of objects. 

In Section 2, we describe the scope of our numerical 
experiments. The model parameters are summarised in 
Section 3. In Section 4, we introduce the simulations and 
present our results in Section 5. In Section 6, we discuss 
these results in light of new observations before giving 
the conclusions in Section 7. 

2. STATEMENT OE THE PROBLEM 
2.1. What is a first galaxy? 

Our interest in galaxy systems with the lowest possible 
mass arises from a desire to understand the chemistry of 
the first stars (Erebel & Bromm 2012; Erebel et al. 2014, 
K13). Our assumption is that the surviving remnants of 
the first galaxies will have small dark halo masses and 
possibly low baryon mass fractions for the reasons out¬ 
lined in §1. Many of these systems will have merged with 
larger galaxies long ago, but some are likely to be falling 
into larger systems today for the first time. 

So what constitutes a first galaxy? This must re¬ 
tain a long-lived stellar system and importantly a confin¬ 
ing dark matter halo.^ Whether low mass stars formed 
alongside or immediately after the first stars is unknown. 
It is likely that the process was patchy, with the very first 
stars forming in rare density peaks, but with an increas- 

^ We recognise that ‘tidal dwarf galaxies’ form a rare class of 
objects that may not have associated dark matter, but generally a 
galaxy is defined by the existence of a confining dark matter halo, 
at least at the time of formation. 
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ing fraction of the universal volume becoming ionized as 
the first stellar generation took hold. 

Here we make the simple assumption that the uni¬ 
verse was enriched to a threshold metallicity of order 
[Fe/H]=-4 before reionization was complete (e.g. Bromm 
& Yoshida 2011). Some of this gas was accreted onto 
low-mass dark matter halos, and star formation began. 
Our threshold metallicity allows us to consider a ‘normal’ 
initial mass function in order to describe how star for¬ 
mation proceeded at these early redshifts (Bromm et al. 
2001; Schneider et al. 2002). Given the stochasticity of 
star formation at our low masses, there may be a signif¬ 
icant amount of low-mass star formation before the first 
supernova event. Some remnants of first galaxies may 
therefore contain stars formed from gas enriched only by 
Population III stars, thereby providing direct informa¬ 
tion on the first stellar yields. 

However, in the systems relevant to our work (ultra¬ 
faint dwarfs, very metal poor DBAs), some of the de¬ 
tected elements (e.g. Sr) can only have arisen from su¬ 
pernova ejecta. Thus, our minimum mass models have 
the more stringent requirement that they can survive 
one or more supernova events during their lifetime, one- 
shot dwarfs those where the gas never recovers from the 
first SNe. In order to experience more than a single 
star-formation episode, a very small galaxy must retain 
enough gas in the face of the energy outflow from a single 
star-forming episode for a second one to occur, perhaps 
at a much later time. A large burst of energy such as 
a supernova could either heat the galaxy ISM and halt 
star formation or completely evacuate the star-formation 
region of gas. This will cut off star formation until such 
time as the halo re-accretes the blown-out gas, assuming 
this is allowed by the environment. Thus, a ‘first galaxy’ 
is one that must retain sufficient gas to form stars in the 
wake of at least a single supernova event. 

This is not a contrived circumstance and we show be¬ 
low (§6) that there may already be observations to sup¬ 
port such a scenario. It is likely that there are situations 
in which only one or two supernovae occur within a single 
star-forming “burst” (Karlsson et al. 2013; Ritter et al. 
2014). Hurley-Keller et al. (1998) identify several star- 
burst epochs in the Carina dwarf with the first of these 
taking place of order ^10 Gyr ago. At that time, the 
surface density of the star formation rate ( 5 *) within the 
core was roughly 10“^ Mq yr“^ kpc“^ over the period 
of the burst, Tq. This timescale is now well established 
in starburst and dwarf galaxies such that < 10 Myr 
(Tolstoy et al. 2009; Sharp & Bland-Hawthorn 2010). 

2.2. Star formation and supernovae 

To test the plausibility of a single supernova event oc¬ 
curring in our models, we consider a simple model of star 
formation. We consider dark-matter potentials with an 
Einasto profile over the virial mass range log[Mvir(M 0 )] 
= 5.5, 6.0, 6.5 and 7.0; we refer to these models as M55, 
M60, M65 and M70 respectively (see Fig. 1 ). If we inte¬ 
grate over the Einasto profile, the total masses for these 
models in dark matter are log[Mtot(M 0 )] = 5.9, 6.4, 7.0 
and 7.5. (These total masses do not include the small 
amount of baryons.) In Fig. 2, we show how the low 
mass halos are expected to grow in mass and radius with 
cosmic time. Halo growth is not treated in our simu¬ 
lations because the time frame of our analysis is short. 



Fig. 1.— Einasto dark matter halos considered here showing 
the relative sizes of the virial and scale radii ai z = 10. The virial 
masses of the four haloes are log[Mvir(MQ)] = 5.5, 6.0, 6.5 and 7.0; 
we refer to these models as M55, M60, M65 and M70 respectively. 
For all halos, the virial mass and radius continue to grow to 2 : = 0 
(see Fig. 2). 


Redshift (z) 



Time (Gyr) 

Fig. 2. — The predicted growth in the virial mass (top) and radius 
(bottom) of the M70 halo at 2 = 10 (dotted line) to the present day 
(see Fig. 1). The shaded band is the redshift range of the damped 
Lyo: systems discussed in later sections. This plot was derived from 
5000 runs of the tree merger code of Parkinson et al. (2008): the 
dashed lines encompass 67% of the predicted halos at each epoch. 

Our results at a given halo mass at z = 10 eorrespond to 
more massive systems at later times. 

Here we show complete results for M65 and M70, and 
partial results for M55 and M60 because these models fail 
to retain sufficient baryons. The initial baryon fraction 
is assumed to be roughly 10 % within the virial radius 
(e.g. Krumholz & Dekel 2012 ). We assume that star 
formation occurs where the gas is densest, i.e. within 
a radius enclosing 50% of the baryons (ro. 5 ^) inside one 
scale radius r^, the scale radius of the Einasto potential 






Lowest Mass Galaxies 


5 



Log[T(K)] 


Log[T(K)] 


Fig. 3. — The variable cooling function, using compositions encountered in the simulations. In each panel, the heavy line is the total 
cooling, the thin red line is cooling due to hydrogen, the green dash-dot line is due to helium, and the metals combined produce the blue 
dotted line. In the single supernova model, the internal metal ratios are fixed by the supernova composition, and are grouped into a single 
function. 




Fig. 4.— A: Dark matter and initial gas density (neutral gas, fi = 1.21, jS = 2.0) profiles for each of the model halos. B: The potential 
functions (—^) and enclosed dark matter masses. Both panels: The dashed, dash-dot, heavy and thin curves are for the M55, M60, M65 
and M70 models respectively. The vertical grey bands cover the radial extent of four scale radii {vs) and four virial radii (rvir), marked 
with four vertical lines, for the four mass models using the same line type for each model. 
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Fig. 5.— A 2D representation of the 3D computational domain. 
Each square represents a cubic domain enclosing the interior re¬ 
gions. The main level of interest is LI, where the mass and energy 
transport during the explosions are measured. A third inner level, 
L2, is used to resolve early epoch (t < 50 kyr) thin shell phases of 
the SNR. 

(see Section 3.1). Over this volume, we find that the gas 
density is within a factor of 6 of the maximum central 
value. 

For the models (M55, M60, M65, M70), the enclosed 
gas mass within the half gas-mass radius (<rs/2) are 
(2.8 X 10^ 9.5 X 10^ 2.5 x 10^, 7.5 x 10^) M© and thus the 
average gas surface densities are (2.1, 2.6, 3.3, 4.2) M© 
pc“^. In the clumpy gas models, there are order of mag¬ 
nitude variations in local gas density. From direct obser¬ 
vation, we know that star formation takes place at these 
moderately low column densities over a range of metal- 
licities down to [Fe/H] ^ -1 (e.g. Bigiel et al. 2008). As 
these authors show, the star formation efficiency is ob¬ 
served to vary wildly below 9M© pc“^ with a median 
value of only a few percent. Little is known about how 
the star formation rate depends on metallicity, but it 
seems likely that it declines with [Fe/H] due to the diffi¬ 
culty of forming H 2 . We stress however that the inferred 
value of 5 * in Carina is comparable to what is observed in 
the outer metal poor regions of spiral galaxies (e.g. Bigiel 
et al. 2008), and this dwarf galaxy somehow managed to 
form long lived stars at much lower metallicities. 

For a star formation efficiency of 2% (Kennicutt 1989), 
the expected total stellar mass over the lifetime of the sys¬ 
tem is (56, 190, 670, 2400) M©. If we assume a Kroupa 
mass function, roughly one fifth of the mass in stars will 
end their lives as supernovae. These estimates are sus¬ 
ceptible to a lot of scatter in the low mass/surface density 
limit (Bigiel et al. 2008). If we scale the projected star 
formation rate per unit area to the largest burst in Ca¬ 
rina (^ 10“^ M© yr“^ kpc“^), the timescales required to 
form these stellar masses are (44, 54, 71, 86) Myr, for a 
half gas-mass radius of (20, 34, 55, 94) pc (see §3.1 for 
further details on the halo properties). If a typical burst 
lasts 10 Myr as in Carina (Tolstoy et al. 2009), we would 
expect multiple bursts over the star-forming lifetime of 
the galaxy. 


So can we consider an isolated supernova event to be 
applicable to all of our models? It is conceivable that 
more than one event can occur within the timescale of a 
single burst Tq. The results given in Table 1 show that 
multiple SN events are rare in the M55 and M60 models, 
are relatively unusual in the M65 model although occur 
half the time in the M70 model. Multiple supernovae 
in more gas-rich M70 models happen half the time, but 
in practice we can consider most of these eventualities 
to be “single supernova events” because one of these is 
highly likely to be a low mass event and physically sep¬ 
arated from the other event. These events are likely to 
enrich their own immediate environs, although substan¬ 
tial mixing can occur between the events at late times 
(Ritter et al. 2014). While different choices of the size 
of the star-forming region or baryon fraction may alter 
these results, we would not expect this to alter the over¬ 
all conclusion of a single supernova event. We note that 
simulations of more massive halos will need to consider 
multiple events within the same burst. 


TABLE 1 

Likelihood of supernova events scaled to Carina star 

FORMATION RATE FOR A TYPICAL BURST DURATION LASTING 10 

Myr; tsN is the median length of time between supernovae 
IN Myr. 



M55 

M60 

M65 

M70 

0 

93% 

80% 

45% 

19% 

1 

6.6% 

18% 

36% 

32% 

>1 

0.2% 

2.0% 

19% 

50% 

tsN (Myr) 

100 

30 

8 
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2.3. The pre-supernova phase 

Massive stars evolve rapidly in the lead-up to the ex¬ 
plosion during which time the stellar wind and radiative 
luminosity can be substantial. This pre-supernova phase 
can be an important factor in removing gas, particularly 
for halos with low baryon fractions. In the worst case sce¬ 
nario, a star at the halo centre can ionize an amount of 
gas, Mion ^ S^mp/ann M© where a is the total recom¬ 
bination coefficient, rup is the mass of the proton, and 
nn is the hydrogen density. Using the Schaerer (2003) 
models at a metallicity of [Fe/H] = -4.0, the mean ioniz¬ 
ing rate averaged over the Kroupa initial mass function 
is S'* ^ 4 X 10^^ phot s“^ for the M55 model, 1 x 10^^ 
phot s“^ for M60, and 3 x 10^^ phot s“^ for M65. Taking 
the worst case scenario where the most massive star falls 
near the centre, the total ionized masses are 1.7 x 10^ M© 
for M55, 6.1 x 10^ M© for M60, 2.2 x 10^ M© for M65 and 
5.5 X 10^ M© for M70. In each case, a large fraction of 
the neutral gas can be ionised. This is likely to have 
a significant impact in reducing subsequent star forma¬ 
tion. But the more likely situation is not expected to be 
so severe. Because the total amount of star formation 
is small in our models, stochastic effects can dominate. 
The median radiation levels expected are more biassed to 
B stars rather than rare O stars, and the ionized masses 
can be much less than 10%. We include the effects of 
pre-ionization in our models below. 
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3. MODEL PARAMETERS 

3.1. Numerical models 

In our experiments, we look at a series of spherically 
symmetric gas distributions in hydrostatic equilibrium, 
detonate a single supernova, and track the ability of the 
potential to retain significant amounts of gas over 25 Myr 
or more. Eor all of our models, we consider the period 
prior to the onset of cosmological reionization when stars 
were forming for the first time. At such early times, sub¬ 
stantial amounts of gas have already begun to settle into 
dark matter halos up to and exceeding Merit (Bullock 
et al. 2001; Ricotti & Gnedin 2005). This is consistent 
with the latest results from deep photometry for ultra 
faint dwarfs where star formation appears to have been 
quenched by reionization in almost all cases (Brown et al. 
2012, 2014; Webster et al. 2015a). 

We consider both clumpy and smooth gas distributions 
within the low-mass halos. Eor the clumpy models, we 
make the important assumption that the confined gas 
has a fractal distribution before the supernova explodes. 
There are two reasons for this assumption. Eirst, the in- 
falling gas becomes shocked and highly inhomogeneous 
close to the dark halo core (Abel et al. 2000). Secondly, 
a high mass star that is destined to become a supernova 
injects mechanical and ionizing energy into the surround¬ 
ing gas prior to the explosion. This energy input drives 
turbulence in the surrounding gas (Elmegreen & Scalo 
2004). 

In looking for possible key factors in the retention, we 
consider the following scenarios: 

• Abundances & cooling: We use chemical abun¬ 
dances that approximate primordial conditions 
starting at [Ee/H] = —4.0. This starting metal- 
licity is conservative: Wise & Abel (2008) show 
that the initial metallicities of the first galaxies 
could have been as elevated as [Ee/H] ^ -3.0 or 
higher. We compute models with cooling, and 
models (7 = |) where we have artificially turned 
off microphysics affecting energy loss (e.g. cooling 
through metal and recombination lines). The cool¬ 
ing function is presented in Eig. 3; the abundances 
are given in Appendix A. 

• Clumpiness: We consider models with quasi¬ 
fractal cool star-forming gas, and comparison mod¬ 
els where the gas is homogeneous and smooth. The 
key point here is that the degree of coupling be¬ 
tween the supernova energy and the surrounding 
medium will depend on whether the medium is 
smooth or structured. 

• Central and off-centre explosions: We consider a 
single 10^^ erg supernova explosion (SNe), located 
either at the geometrical centre of the halo poten¬ 
tial, or off-centre at a radius that encloses half the 
mass contained inside one scale radius of the poten¬ 
tial. Supernovae that occur away from the centre 
have an asymmetric environment. 

This gives, in combination, 8 scenarios for each of the 
halo masses considered. In each case, the mass and en¬ 
ergy contained in the inner scale radius is integrated at 


25 kyr intervals and plotted over the course of the simu¬ 
lation to obtain estimates for the masses retained at late 
times. 

In the Fyris Alpha models presented here, we use a 
fixed dark matter potential that dominates everywhere. 
Euture models will include self-gravity in the gas. There¬ 
fore we limit the initial baryonic mass inside a scale ra¬ 
dius (rg) to be 10% of the dark matter mass inside that 
radius. This results in an overall baryonic fraction out to 
the virial radius (rvir) of /b ~ 12%, less than the canon¬ 
ical /crit ~ 17% global baryonic mass fraction (Planck 
Gollaboration et al. 2014). Our choice of /b is a conser¬ 
vative assumption because a lower gas budget generally 
assists the supernova in accelerating and heating the gas. 
If the models retain gas with these lower mass fractions, 
they can be reasonably expected to do so with more gas 
in place. 

3.2. Dark matter halo 

In recent years, high resolution GDM models have iden¬ 
tified the need for a halo density profile that has an addi¬ 
tional parameter beyond those used in the NEW profile 
(Springel et al. 2008; Hayashi et al. 2004). The profiles 
seen in the latest simulations display local density deriva¬ 
tives that are power laws of radius, unlike the NEW pro¬ 
file that has asymptotic behaviour at small radius. They 
invoke a generalised exponential function that has this 
property, and note that it has been used in the past by 
Einasto (1965) in the context of galactic halo spatial den- 
ities^. 

We use the Einasto function, a generalised exponential, 
to form a dimensionless potential defined by the parame¬ 
ters a, a halo concentration factor Xvir (often referred to 
as c) and the density contrast of the halo A (^178) with 
respect to the universal background (Nichols & Bland- 
Hawthorn 2009). We stress that these parameters must 
be specified at a fixed redshift z because of their strong 
cosmic evolution (Power et al. 2003; Duffy et al. 2008). 

We consider four halo mass models (M55, M60, M65, 
M70) defined in Section 2.2. The total halo masses are 
(3.0, 2.7, 2.5, 2.3) times the virial mass respectively (e.g. 
Fig. 4). The halo and gas properties are summarised in 
Table 3. Our normalization uses a density function value 
of 1/e at the scale radius, r^. This gives a slightly dif¬ 
ferent set of scaling constants, but recovers exactly the 
same physical quantities as the |-scaling used in Nichols 
& Bland-Hawthorn (2009). This change was motivated 
by wanting to compare with dimensionless isothermal po¬ 
tentials used by Sutherland & Bicknell (2007), which also 
obtained a scaled density value of 1/e at the scale radius. 
(We refrain from using the popular ‘core’ radius because, 
unlike the isothermal potentials, the generalized Einasto 
functions do not have a distinct core region; here we use 
rs throughout.) 

Our adopted halo parameters are based on those of 
Duffy et al. (2008) and Springel et al. (2008). As their 
fits to the simulated halos are measured over a relatively 
low redshift range, we extrapolate these properties back 
to a redshift of z = 10. At the time of writing, we became 
aware of new work (Gorrea et al. 2015) where the param- 

^ This has the same functional form as the projected Sersic pro¬ 
file but the 2D function does not transform to the same underlying 
radial density distribution. 
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TABLE 2 

Model Parameters: the cosmology is fixed by the halo formation redshift, z, the normalized Hubble constant h = Ho/100 

KM S“^ MpC“^, the dark ENERGY DENSITY AND THE DARK MATTER DENSITY Qm- THE EINASTO HALO PROPERTIES (SEE TEXT) ARE 
FIXED BY THE CONCENTRATION PARAMETER, O, AND THE HALO DARK MATTER DENSITY CONTRAST WITH THE AMBIENT MEDIUM, A, DEFINED 
AT REDSHIFT Z. THIS DENSITY CONTRAST CAN BE EXPRESSED AS Ac AVERAGED OVER THE VIRIAL RADIUS AND NORMALIZED TO THE MEAN 


CRITICAL DENSITY pc, OR AS A^^ 

NORMALIZED TO THE 

UNIVERSAL MASS DENSITY (pu)- 

The halo potential is defined by its scale 

RADIUS Ts, THE MEAN DENSITY AND MASS 

WITHIN THIS 

RADIUS Ps AND Mg 

. The halo VIRIAL mass is Mvir WITHIN THE VIRIAL RADIUS 

^vir = M 300 IS THE MASS 

WITHIN 300 PC AND Mtot IS THE MASS INTEGRATED TO 

1 INFINITY. For the baryons, Cs and Cs^iso ARE 

THE ADIABATIC (PHYSICAL) AND ISOTHERMAL SOUND 

SPEEDS IN THE GAS 

WITH MEAN TEMPERATURE Tgas- THIS GAS HAS A CENTRAL 

DENSITY (number DENSITY) OF 

Pgas,0 ('Rif,o) WITH A 

VIRIAL (total) mass OF Mgas,vir (Mgas,tot)- NOTE THAT THE INITIAL BARYON 

FRACTION fb ^ Mgas,s/Ms 

~ ^gas,vir/^vir ~ 10% AND THE RADIUS 

1 CONTAINING 

HALF THE INITIAL GAS WITHIN Tg IS ro.5^. 




Cosmology 





z 

h 

Ha 

Hm 



10.0 

0.7 

0.7 

0.3 


a 

Ac 

Pc(g/cm3) 

Ati 

{PM>(g/cm3) 


0.18 

177.5 

3.68E-027 

591.7 

l.llE-027 




Dark Matter 





Ts 

Ps 

a^vir 

■^vir 


Model 

(pc) 

(g/cm^) 


(pc) 


M55 

33.2 

1.02E-023 

5.98 

198.5 


M60 

55.0 

7.73E-024 

5.29 

291.3 


M65 

91.3 

5.87E-024 

4.69 

427.6 


M70 

151.3 

4.48E-024 

4.15 

627.7 



Ms 

Mvir 

M 300 

A^tot 


Model 

(Mo) 

(Mo) 

(Mo) 

(M©) 


M55 

5.64E+04 

3.16E+05 

4.10E+05 

7.90E+05 


M60 

1.94E+05 

l.OOE+06 

1.02E+06 

2.73E+06 


M65 

6.73E+05 

3.16E+06 

2.36E+06 

9.44E+06 


M70 

2.34E+06 

l.OOE+07 

5.08E+06 

3.29E+07 




Baryons 





Cs 

Cs, iso 

Tgas 

A^gas,vir 


Model 

(km/s) 

(km/s) 

(K) 

(M©) 


M55 

2.70 

1.91 

5.33E+02 

3.87E+04 


M60 

3.90 

2.76 

l.llE+03 

1.19E+05 


M65 

5.63 

3.98 

2.32E+03 

3.66E+05 


M70 

8.16 

5.77 

4.86E+03 

1.15E+06 



ro.5g 

Pgas,0 

^H,0 

ATgas© 


Model 

(pc) 

(g/cm®) 

(cm-3) 

(M©) 


M55 

20.8 

4.09E-23 

17.2 

5.63E+03 


M60 

34.5 

3.11E-23 

13.1 

1.94E+04 


M65 

57.1 

2.36E-23 

9.91 

6.73E+04 


M70 

94.8 

1.79E-23 

7.55 

2.34E+05 


eters are extended to back to z = 20 but these appear to 
be consistent with our estimated values, particularly for 
the important M65 and M70 models. 

3.3. The cold interstellar medium 

In our simulations, we introduce a non-uniform 
medium to study the influence of inhomogeneity on the 
dynamical interaction between the supernova ejecta and 
the surrounding gas. To assess the importance of in¬ 
homogeneity, we run both homogeneous (smooth) and 
inhomogeneous (clumpy) models and compare the out¬ 
comes. Following Sutherland & Bicknell (2007), to es¬ 
tablish a non-uniform medium, we make use of an anal¬ 
ogy with a turbulent medium. Rather than attempt to 


review this huge topic, we refer the reader to astrophysi- 
cally oriented reviews as a starting point for background 
material (Elmegreen & Scalo 2004; Scalo & Elmegreen 
2004). 

We represent the non-uniform properties of the turbu¬ 
lent medium using three statistical characteristics: the 
variance, the intermittency (by using log-normal dis¬ 
tributions), and a self-similar power-law structure. We 
achieve both the property of intermittency and power- 
law structure simultaneously with an iterative scheme 
first developed by Lewis & Austin (2002); see Appendix 
B for more details. Hence the initial distribution of the 
ISM that we employ should be regarded as a physically 
motivated generalization of a homogeneous model. 



Lowest Mass Galaxies 


9 


We further simplify our turbulent medium by neglect¬ 
ing the velocity structure, and focus on choosing sta¬ 
tistical parameters to describe the density alone. This 
is justified numerically by the relatively small turbulent 
velocities expected when compared to the very large ve¬ 
locities found in our global SNe-ISM interaction. The 
velocities observed in typical warm ISM conditions fall in 
a range of transonic to mildly supersonic values, Mach 1- 
5 (e.g. Heiles 2004). At temperatures at or below 10^ K, 
this corresponds to velocities ^ 4 km s“^. In our 0.1- 
1 kpc simulations over 10^ yr timescales, the resulting 
displacements amount to only a few cells at most, and 
are insignificant compared to those of the energy bubble 
generated by the supernova outburst, where velocities of 
many hundreds of km s“^ occur. Consequently, we do 
not impose a turbulent velocity field on the cold ISM 
medium. 

3.4. Hydrodynamic calculations 

The modelling was carried out using the Fyris Alpha 
code (Sutherland 2010; Bland-Hawthorn et al. 2007). It 
solves ideal Euler hydrodynamic flows of gas with cool¬ 
ing and an equation of state appropriate to astrophysical 
gasses, in the presence of a fixed gravitational field. The 
code uses a third order semi-Lagrangian shock-capturing 
method that is robust and well suited to the extreme 
temperatures and Mach numbers encountered here. 

The computational domain consists of a three dimen¬ 
sional nested series of cubic meshes, as illustrated by the 
schematic in Fig. 5. Each level is a factor of 3 smaller 
in physical domain and higher resolution than its outer 
containing level, covering the central third of the enclos¬ 
ing level. Each level is 216^ cells with an overall effective 
resolution of 1944^ cells. 

A relatively low resolution outer level, LO, covers the 
halo out to ±6 scale radii, and provides an outer reservoir 
and boundary condition for the main level of interest, 
LI. The outer reservoir reduces or avoids completely 
any boundary artefacts on LI due to gravity bringing 
new gas onto the main level when infall conditions ap¬ 
ply, keeping level LI as physically consistent as possible 
at all times. LO contains the virial radius. Velocities 
at the outer boundary are extremely small at all times 
considered here and the outer LO boundary is essentially 
fixed. 

The main level, LI, covers a region of just over ±2 scale 
radii, and is the primary region where the movement and 
evolution of mass and energy are measured. Early epochs 
require higher resolution to resolve thin spherical shell 
radiative blast-waves, so a third, inner injection level, 
L2, is added for to resolve the early phases of the SNe. 
After 25 kyr, L2 is disabled and the remainder of the 
simulations are carried out on LO and LI. 

For each halo mass, the scale radius varies and so the 
resolution is variable. The highest resolution, on level 
L2, is Ax ~ 0.25, 0.37, and 0.62 pc per cell for the M55, 
M60 and M65 models. Thus the effective mass resolution 
is ^ 10^ Mq / 216^ or of order 10“^ Mq; in practice, cell 
densities can vary by more than 4 orders of magnitude 
across the grid. 

The remaining key features of the calculations can be 
summarised: 

• A fixed gravity potential, interpolated locally by 


3rd order PPM interpolation to match the inter¬ 
polation of the hydro dynamical algorithm, taking 
conservative potential differences to compute mean 
accelerations. 

• SNe simulated by depositing 10^^ erg of energy as 
internal gas energy in a smoothed spherical region, 
as small as possible while retaining spherical ge¬ 
ometry well enough to prevent gross non-spherical 
‘pixel’ errors, typically 6 cell radius on L2, corre¬ 
sponding to 1.5—3.7 pc. This then takes the form of 
a pulse of low density hot over-pressure gas which 
then expands, converting internal energy to kinetic 
energy and doing work against the potential. En¬ 
ergy is also lost via cooling. 

• To allow for significant changes in the local metal- 
licity in the models as the metal rich-ejecta prop¬ 
agate through the models, a generalised cooling 
model has been adopted that allows the local com¬ 
position of each cell to be determined during the 
simulation. Instead of a single cooling function and 
uniform molecular weight, the cooling was sepa¬ 
rated into three components, representing hydro¬ 
gen, helium and heavier metals, and ionization was 
taken into account. 

The metal-rich supernova ejecta were followed using 
a scalar tracer variable, advected passively, representing 
the supernova material, and hence the local helium frac¬ 
tion and metal fraction. The local enrichment of the 
gas was thus computed at each time step. Furthermore, 
the ionization fraction of the components, precomputed 
with the MAPPINGS IV code (Allen et al. 2008; Suther¬ 
land et al. 1993) as a function of temperature, allowed 
the equation of state to vary with the mean molecular 
weight, ensuring more accurate shock temperatures, and 
allowing for a wider range of cooling timescales caused by 
strongly varying composition. Fig. 3 shows the resulting 
computed cooling functions for a range of encountered 
compositions on the grids. 

We acknowledge that there are large uncertainties in 
the yields from essentially zero metallicity stellar super¬ 
nova models. The yield of metals, particularly oxygen, 
could have an impact on the outcome, lower mass super¬ 
novae giving less enrichment and cooling, while higher 
mass stars could give more cooling. More detailed fu¬ 
ture work incorporating a range of supernova progenitor 
masses will investigate the impact of different supernova 
assumptions. By allowing for the ionization and compo¬ 
sition to vary from cell to cell, we obtain more accurate 
species columns, for example, for neutral and ionised hy¬ 
drogen. 

Our focus is on the hotter atomic and ionised phases 
of high pressure gas, and the cooling therein, as that gas 
is volume filling and will have the primary effect on the 
global dynamic evolution if any. Gas in the molecular 
phase is likely to be confined to the centres of the dens¬ 
est cores, if present at all. As the outflows are all warm 
(> lO^K), they are expected to be atomic or ionised in 
any case. In the face of uncertainty in the formation of 
H 2 under primordial conditions and little or no informa¬ 
tion on GO cooling, we leave molecular cooling out. It 
would be essential to include it if we were tracking the 
(much longer) time required for the hot gas to re-cool 
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back to star-forming temperatures ^ 100 K, but here 
we are focussed on the initial halo gas retention problem. 
We defer the star formation problem to more substantial 
investigations. 

That we neglect molecular cooling in these simulations 
can be considered a conservative approximation, but we 
expect it to have little effect on the global dynamics of 
gas retention. More than 99.9% of the internal energy 
is lost by the time the gas cools back below 10^ K from 
early temperatures of more than 10^ K, so only a tiny 
fraction of the SNe energy is available to the molecular 
cooling in the present simulations. 

In Table 4, we calculated the thermalization efficiency 
for a subset of the models as a check that our simulations 
correctly handle cooling. The thermalization efficiency 
was estimated by finding the ratio of the total energy at 
the time when gas is first leaving the grid to the total ini¬ 
tial energy. Our efficiencies should be considered upper 
bounds as it is likely that some cooling will occur at later 
times. From Thornton et al. (1998), we expect the ther¬ 
malization efficiency to be ^ 10%, which is consistent 
with the results of our models. 

3.5. Pre-ionization from the SN progenitor 

The pre-ionization of the gas by the supernova precur¬ 
sor plays a key role in the loss or survival of gas after 
the supernova explosion, particularly in low-mass dark 
matter halos. We consider in the detail the impact of an 
ionised HII region that forms around the progenitor star 
prior to the explosion. The key parameters of this model 
are the stellar properties (Teff, L) and the associated stel¬ 
lar atmospheres. These give the ionising photon flux as 
a function of time up to the SN event. Transforming the 
grid coordinates to spherical coordinates around the star, 
radial summations to approximate the optical depth inte¬ 
grals are performed, resulting in thermal and ionization 
structures calculated from the MAPPINGS /Fionization 
code. These are then mapped back to cartesian coordi¬ 
nates and used in the Fyris fluid computations. It was 
necessary to perform this remapping for every timestep 
in the fluid dynamics simulation in order to track the 
evolution in opacity. This slowed the code by a factor of 
5-10. We considered both clumpy and smooth external 
media with and without gas cooling. 

We model the supernova progenitor with a low metal- 
licity star of 25 Mq. Higher mass stars produce stronger 
ionising intensity and therefore assist gas loss by over¬ 
heating the gas prior to the supernova explosion. But 
since 80% of all supernova events in a conventional ini¬ 
tial mass function have lower mass and therefore weaker 
ionising radiation fields, we consider our choice of the up¬ 
per mass limit to be conservative. We use the evolution¬ 
ary tracks of Meynet & Maeder (2002, hereafter MM02) 
for a rapidly rotating, low-metallicity massive star. The 
25 Mq track was missing from that work, and so we have 
interpolated the main sequence phase of this model from 
the other stellar models. In Fig. 6A, we show the inter¬ 
polated values of the luminosity and temperature for the 
25 Mq model at [Fe/H] ^ -3.7 (MM02). These are the 
closest published tracks for our required starting metal- 
licity of [Fe/H] = -4.0. We used the tabulated MS life¬ 
time of 6.1 Myr; the interpolated temperature-luminosity 
track is shown in Fig. 6B. 

In Fig. 7, we show the impact of the total ionising 


budget during the pre-supernova phase for our adopted 
stellar model and a perfect black body of the same sur¬ 
face temperature. The panels A, C and E refer to the 
models with the progenitor at the centre; panels B, D, F 
are for the off-centred models. This gives some indication 
of the sensitivity to the stellar models. If the thin lines 
exceed the thick lines, then a Stromgren sphere is possi¬ 
ble within the radius, and the stellar radiation would be 
incapable of extending beyond that, which is the initial 
case for all halos. However as the gas responds to heating 
and evolves, the density is smoothed and falls on average 
due to outflows. The quadratic dependence of recom¬ 
binations means that the available recombinations falls 
rapidly, and once the photon curves exceed the recombi¬ 
nation curves, the photoionised zones can grow rapidly 
and ionise more of the halo. 

In the upper panels, the M60 halo will be easily fully 
ionised after the first 2 Myr. The middle panels show the 
M65 models are marginal if the photons are uniformly 
spread throughout the halo, and capable of ionising the 
inner core region and possibly beyond before the super¬ 
nova occurs. The actual result will depend on the details 
of the radiative transfer, whereas this is an indicative 
integral view of the processes which does not take the 
geometry of the ionised volume into account. The lower 
panels for M70 indicate that the stars are insufficient to 
photoionise a significant volume of the halo. 

These indications are broadly supported by the simu¬ 
lations: M60 is rapidly ionised, M65 is able to break out 
of the core radius before the SNe, and the ionised region 
is contained in the M70 models. 

Stellar winds. To avoid the complexities of rapid post- 
main sequence evolution, we assume that the star simply 
explodes as a supernova at the end of the main sequence. 
Ekstrom et al. (2006) showed that mass loss during the 
main sequence in similar stars at very low metallicity may 
only amount to 1% of the initial mass, and Kudritzki 
(2002, 2005) showed that the stellar wind luminosities of 
low metallicity ([Ee/H]< —4.0) O stars are typically 10^^ 
erg s“^ or less, several magnitudes below an equivalent 
solar metallicity star. We combine these values to get a 
mean wind velocity, for a main sequence mass loss of 1%: 


1034erg s“ 


^MS 

6.1Myr 


This wind is included in the pre-supernova evolution, but 
the small mass flux and low ram pressure meant that it 
had negligible effect on the simulations. Once the HII re¬ 
gion pressurised the surrounding region, the stellar wind 
was unable to blow a wind bubble of any significant size. 

The key parameter for these near zero metallicity O 
stars is their unusually high effective temperatures for 
a given mass. At 25 Mq, we estimate the initial tem¬ 


perature to be 48,000 K , compared to a similar mass 
star at solar metallicity of only 40,000 K (see for exam¬ 
ple Meynet & Maeder (2000) and more recently (Georgy 
et al. 2012)). This means the star will be more effective 
as a source of ionising photons than the equivalent star 


in the solar neighbourhood. To convert the evolutionary 
tracks, Teff(t) and T(t), to an ionising photon luminos¬ 
ity (see Eig. 8), we used the ATLAS9 atmospheric grid 
(Castelli & Kurucz 2004) which uniquely contains a set of 
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Fig. 6. — A: Interpolated luminosity and temperature for the 25 Mq model in the MM02 grid. The zero age main sequence is labelled 
ZAMS and the end of the main sequence as MS. B: Evolutionary tracks for the four most massive rotating (300 km s“^) star models in 
MM02 at [Fe/H] -4. We adopt a 25 Mq star model interpolated from these tracks with cubic splines using the endpoints in A. For this 
star, the time span from the zero age main sequence (ZAMS) to the explosion is 6.1 Myr. 


atmospheres for [Fe/H]= —4.0 with alpha element abun¬ 
dance enhancements, perfectly matching the initial gas 
composition in our simulations. The interpolated stellar 
temperature, luminosity and ionising fluxes are shown in 
Fig. 8. The derived ionizing spectrum is shown in Fig. 9; 
for comparison, a photon flux assuming a blackbody in¬ 
stead of the stellar atmosphere is included. 

4. SIMULATIONS 

In the absence of considering the full initial value prob¬ 
lem (i.e. infall, self gravity), we assume that sufficient gas 
has already settled onto a spherical dark-matter halo, 
and that it has dissipated energy and come to rough dy¬ 
namical equilibrium. Fundamentally, an isothermal hy¬ 
drostatic equilibrium in a gravitational potential, 4>, is 
only possible with a uniform medium, with the require¬ 
ment that the local pressure gradients oppose the local 
potential gradient, that is: 

14 ^ = /ON 

p dr dr ^ 

where P and p are the gas pressure and density at a 
radius r. As the potential is smooth, the pressure gradi¬ 
ents must be well behaved. When the density is smooth 
also, a solution for a given temperature is possible, with 
the temperature defining the pressure scale height. How¬ 
ever, when the density is no longer locally smooth, cor¬ 
responding local temperature variations are required to 
give the well-behaved pressure gradients needed for hy¬ 
drostatic equilibrium. Without a single global tempera¬ 
ture, a range of scale heights exist over the whole domain, 
and a single structure scale - and hence equilibrium - 
is not strictly possible. However, in practice it is close 
to equilibrium, and the timescales for change are much 
longer than those that describe the star formation in our 
models. 

When gas falls into a gravitating potential well, it un¬ 
dergoes complex motions as it settles towards dynami¬ 


cal equilibrium. Some of the kinetic and potential en¬ 
ergy is likely to be converted to thermal energy which 
is then lost through radiative cooling. If there were 
no losses, after adiabatically virialising in the halo po¬ 
tential, the gaseous medium would have an initial tem¬ 
perature, Tinit, related to the virial temperature of the 
halo such that Tinit ^ ^vir- For the adiabatic case, Tvir 
= for which k is Boltzmann’s constant, 

p is the mean molecular weight of the gas and mu is the 
mass of the H atom. This is before any star formation 
takes place. 

We now introduce the common (3 notation used in 
the presence of complex multiphase gases (Cavaliere & 
Fusco-Femiano 1976; Arnaud 2009). An initial (quasi) 
hydrostatic gas distribution depends on an effective tem¬ 
perature Teff = Tgas + Tturb ^ Finit arising from a com¬ 
bination of thermal and turbulent support (plus other 
terms due to magnetic fields which we neglect at present). 
With radiative losses, the initial thermal pressure is lost, 
and we are left with turbulent support only. Thus, if 
we assume initial equipartitition between thermal and 
turbulent pressure, Teff = Tturb ~ 0.5Tvir such that 

P = Tvir/Teff = 2.0. 

For all of our models, we adopt a /d = 2.0 gas dis¬ 
tribution to recognize that the baryons must undergo 
some cooling before the first stars have formed. It is 
conceivable that P is larger but such an assumption 
would concentrate the gas closer to the centre of the 
potential, thereby assisting gas retention. Indeed, Wise 
et al. (2008) find that baryonic matter is more con¬ 
centrated than the dark matter in their evolving halos 
(Mtot= 10^ Mq). These authors use an adaptive mesh 
code to track the development of supersonic turbulence 
driven by accretion over many orders of magnitude in 
physical scale. In combination with line cooling, the su¬ 
personic turbulence leads to central densities that are 
several orders of magnitude higher than what we can 
achieve with Fyris. 
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Fig. 7. — Comparison of the instantaneous photon rate (heavy lines) compared to the total available recombinations as the models 
progress in the pre-supernova phase. The lefthand figures (A, C, E) correspond to the centred models; the righthand figures (B, D, F) refer 
to the off-centred models. The two heavy lines represent a pure black body model (BB) and an extreme low metallicity stellar atmosphere 
(Star), and the difference between them gives an indication of the magnitude of any uncertainties in the atmosphere models. Available 
recombinations are computed by integrating neuJa(T 4 ) in cm“^ s“^ over the respective volumes where ol is the recombination coefficient 
summed over all possible transitions and T 4 is the gas temperature in units of lO^K. This yields the total number of recombinations per 
second if all of the gas can be ionised, ignoring optical depth and geometry of the ionised volume. 
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Fig. 8 .— Interpolated L{t) [long dash], Teff(t)[short dash] and 
ionising photon fluxes, [heavy lines] for M25. Fstar assumes stellar 
atmosphere models, and a comparison is shown, Fbb? which is the 
flux assuming a simple black body approximation to the stellar 
spectrum. The star emits between 5 x 10^® and 7 x 10^® ionizing 
phot s“^ over its lifetime. 



Fig. 9. — Stellar spectrum (solid line) of our adopted 25 Mq star 
model using an ATLAS9 ZAMS stellar atmosphere at [Fe/H] -4 
(Castelli & Kurucz 2004). The dashed line is the blackbody curve 
for the same surface temperature Tgff = 48 000 K. 


While we have the option to modify [3 to increase the 
central gas densities, this would only assist our ability 
to retain gas after the supernova. Our adoption of a 
lower constant p value has the benefit of reducing any 
central concentration differences in the models, leaving 
just the global gravity potential as the key mass model 
difference. The resulting temperature distributions for 
the M55-M70 models are discussed below. 

Initial state of the gas. 

We adopt a mean equilibrium distribution for the gas 
density given by p = p(0) exp(<h)^ for which ^ describes 
the dark halo potential (Fig. 4, right panel). This is 
described in detail in Appendix B. The use of f3 allows 
for the gas scale length to be different from Tvir (Fig. 4, 


left panel). For completeness, note that the dynamical 
dispersion of the dark halo is a = Pcg. For a given tem¬ 
perature and composition/ionization state, P (x p. So 
for homogeneous models, there is a uniform tempera¬ 
ture everywhere such that the local density p is equal 
to the mean density {p/{p) = 1.0). For fractal models, 
p/{p) = f{x,y,z) for which / is the normalised fractal, 
giving a distribution T (x Pjp that is the same for all 
the models using the same fractal modulation. The dis¬ 
tribution / has a mean value of unity and a variance of 
5 (see Appendix B). 

The projected 2D gas distribution (density, pressure, 
temperature) for the inhomogeneous M60 model is shown 
in Fig. 10. The temperature histogram for Fig. lOD, and 
for all models, is given in Fig. 11. Across the models, 
the global mean gas temperature spans the range 500- 
2000 K, with a small gas fraction extending below a few 
hundred Kelvin, compatible with a star forming region. 
The bulk of the ISM tops out at about 10^ K and thus 
covers the range of temperatures seen in typical warm 
and cold neutral ISM media (Sternberg et al. 2002). 

To verify the equilibrium of the initial conditions, we 
ran adiabatic test models with smooth gas in hydrostatic 
equilibrium, and the initial distribution remained steady 
for the ^200 Myr, much longer than the 25 Myr timescale 
of our supernova model. An adiabatic clumpy test model 
could not achieve perfect equilibrium, and evolved to¬ 
wards the smooth solution, smearing out the clumps on 
timescales of 

7-cross ^ rs/cs ^ 100 pc/1 km s“^ ^ 100 Myr, (4) 

where is the scale radius of the Einasto halo; Cs is 
the thermal sound speed for a gas temperature of ^ 10^ 
K. This remains a significantly longer timescale than the 
supernova evolution timescale, so we are justified in as¬ 
suming equilibrium during the main evolution models, 
and a differential experiment where we can compare di¬ 
rectly fractal and smooth models that have the same 
mean properties averaged over the distribution is valid. 

Energy injection. To model the energy injection from su¬ 
pernovae, we insert a bubble of hot gas with the equiva¬ 
lent energy of 10^^ erg, on the assumption that the super¬ 
nova energy has been thermalized on a scale much below 
our parsec-scale resolution. This high pressure bubble 
expands and converts internal thermal energy to kinetic 
energy and radiative losses during the expansion. 

For the smooth models with cooling, it is important 
to provide additional smoothing to the boundary of the 
injection region to avoid the grid cell errors that seed the 
well known “carbuncle instability” (cf. Sutherland 2010). 
This is especially strong in the case of smooth radiative 
models because of the formation of radiative thin shells 
which are particularly susceptible to this. We reduce the 
impact of these effects by using an injection region with 
a 3-cell radius smoothed with a gaussian kernel (2 cell 
FWHM). 

We locate our injection region either at the centre of 
the potential or off-centred at the half gas-mass radius 
(Fo. 5 g) listed in Table 3 found from integrating the mean 
density curves in Fig. 4. In Fig. 10, the injection region 
is shown for the MOO off-centred clumpy model. 

We follow the evolution of the SN for Tq = 25 Myr 
taking snapshots at 25 kyr intervals. This timescale is 
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Fig. 10.— The initial distributions of density, pressure and temperature (level LOl) for the MOO off-centred supernova model; the 
distributions are similar for M55, M65 and M70. While the density and temperature structure are fractal, the overall pressure is constant 
at the start. Each panel shows a slice through the middle of the halo and the quantity shown is scaled logarithmically; the circle has a 
radius of r^, the scale radius. A: relative column density over the range -4.0 (white) to 0.0 (black). B: relative density over the range 
-5.0 (white) to 0.0 (black). The off-centred SN is visible as a white dot. C: relative pressure ranging from -5.0 (white) to 0.0 (black). D: 
temperature ranging from 1.0 (white) to > 4.0 (black). 



Fig. 11.— The initial gas temperature distribution for the mate¬ 
rial inside one scale radius, rs, for the models with fractal clumpy 
media. Voids have higher temperatures and occupy a larger volume 
than the cooler, denser clouds. The material between 10^'^ K and 
10® K arises in the small spherical region that is evacuated prior 
to the supernova injection onto the grid. There the local density 
is reduced to 0.01 of its original values. The fractional volume this 
occupies remains small compared to the entire region inside the 
scale radius. See text for details. 

now well established in starburst and Carina-like dwarf 
galaxies (Sharp & Bland-Hawthorn 2010; Tolstoy et al. 
2009). We analyse our results by summing the dynam¬ 


ically cold and total mass, and internal kinetic energy, 
within one scale radius of the central potential at each 
time step. 

5. RESULTS 

5.1. Pre-ionization phase 

Our goal is to establish the low-mass limit of a spheri¬ 
cal dark-matter halo that can retain baryons after a sin¬ 
gle supernova event. But hrst we consider the impact 
of the progenitor star on the surrounding gas prior to 
the supernova explosion. The mechanical wind energy is 
completely dwarfed by the star’s radiative phase for any 
progenitor of a supernova. 

We consider how the UV propagates in each of the 
four halos for both clumpy and smooth media which 
undergo cooling. This is a complex process and so 
we include access to 3D movies that reveal how the 
ionization front develops. These can be found here: 
http://miocene.anu.edu.au/smahgalaxy. We show a 
time sequence of the full pre-ionization phase and the 
subsequent SN explosive phase in two ways: (i) hxed ref¬ 
erence frame; (ii) rotating reference frame. There are two 
other aids to visualization. First, we include half-tone hg- 
ures that show critical time steps for the different models: 
M60 (Figs. 12, 13), M65 (Figs. 14, 15), and M70 (Figs. 
16, 17). Secondly, we include colour composite hgures in 
Fig. 18 that directly compare the information embodied 
in Figs. 12 to 17, i.e. the neutral and ionised components 
as the nebula develops before and after the SN event. 

A summary of gas mass and energy retention after the 
pre-ionization phase is given in Table 5 and Table 6 re¬ 
spectively. The initial gas masses inside two cardinal 
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TABLE 3 

Naming convention for all numerical simulations 


Models without pre-ionisation 





Central Explosions 


Off-Centre Explosions 




Smooth 

Clumpy 

Smooth 

Clumpy 

Halo Mass 

Adiabatic 

Cooling 

Adiabatic 

Cooling 

Adiabatic Cooling 

Adiabatic 

Cooling 

105.5 

M© 

M55CSA 

M55CSC 

M55CCA 

M55CCC 

M550SA M550SC 

M550CA 

M550CC 

10® 

M© 

M60CSA 

M60CSC 

M60CCA 

M60CCC 

M60OSA M60OSC 

M60OCA 

M60OCC 

106.6 

M© 

M65CSA 

M65CSC 

M65CCA 

M65CCC 

M650SA M650SC 

M650CA 

M650CC 

10'^ M© 








Models with pre-ionisation 




Central Explosions 


Off-Centre Explosions 




Smooth 

Clumpy 

Smooth 

Clumpy 

Halo Mass 

Adiabatic 

Cooling 

Adiabatic 

Cooling 

Adiabatic Cooling 

Adiabatic 

Cooling 

105.5 

M© 








10® 

M© 




M60CCH 



M60OCH 

10®. 5 

M© 




M65CCH 



M650CH 

10^ 

M© 




M70CCH 



M70OCH 


TABLE 4 

Thermalisation efficiency in the M60 and M65 models 



Central Explosions 

Off-Centre Explosions 


Model 

Smooth 

Clumpy 

Smooth 

Clumpy 


Time 

Adiabatic Cooling 

Adiabatic Cooling 

Adiabatic Cooling 

Adiabatic Cooling 

M65 

100% 8.1% 

100% 11% 

100% 9.8% 

100% 

11% 

M60 

100% 11% 

100% 15% 

100% 11% 

100% 

16% 


TABLE 5 

Percentage of initial mass retained inside Ts and rvir for cooling models. The time at 0 Myr is when the SN event takes 

PLACE AT WHICH POINT THE PRE-IONIZATION PHASE ENDS. ENTRIES IN BOLD FALL BELOW 50%. 


Model 

Inside Vg 

Inside rvir 

Time^ 

Central 

Off-Center 

Central 

Off-Center 

M60H 

Mg,s: 1.91E+04 M© 

-^g,vir • 

1.22E+05 M© 

0 Myr 

12.7 

16.1 

97.6 

97.1 

5 Myr 

0.03 

0.06 

89.3 

88.4 

25 Myr 

0.01 

0.21 

0.98 

6.10 

M65H 

Mg,s: 6.61E+O4M0 

Mg,vir: 

3.74E+05 M© 

0 Myr 

50.3 

78.8 

100.0 

99.4 

5 Myr 

2.42 

41.7 

99.7 

97.9 

25 Myr 

13.8 

33.9 

81.5 

73.3 

M70H 

Mg,s: 2.38E+O5M0 

Afg^vir • 

1.16E+06 M© 

OMyr 

102 .6* 

98.2 

100.0 

100.0 

5 Myr 

106.2* 

92.8 

100 .2* 

100 .1* 

25 Myr 

121.5* 

107.0* 

102 .6* 

101 .2* 


Bold: < 50% retained 

L 0 Myr = tsNe? =6.1 Myr from start of simulation 
*: Increased Mass after some infall 
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TABLE 6 

Percentage of initial energy retained inside Ts and rvir for cooling models. The time at 0 Myr is when the SN event takes 

PLACE AT WHICH POINT THE PRE-IONIZATION PHASE ENDS. ENTRIES IN BOLD FALL BELOW 50%. 


Model 

Inside Vg 


Inside rvir 

Time"^ 

Central Off-Center 

Central 

Off-Center 

M60H 

Eg: 4.52E+48 ergs 

Evir: 2.75E+49 ergs 

Pre-SNe 

237.7 

299.1 

1725 

1647 

Post—SNe 

1.74E+04 1.79E+04 

4365 

4374 

5 Myr 

1.90 

1.00 

1242 

1179 

25 Myr 

0.06 

0.31 

1.93 

9.32 

M65H 

Eg: 3.19E+49 

ergs 

Evir: 1.73E+50 ergs 

Pre-SNe 

327.4 

318.1 

325.3 

402.4 

Post-SNe 

2613 

2596 

722.1 

797.3 

5 Myr 

2.71 

50.7 

232.7 

255.4 

25 Myr 

20.1* 

29.6 

101.4 

91.5 

M70H 

Eg: 3.45E+50 

ergs 

Evir: 1.63E+51 ergs 

Pre-SNe 

92.2 

89.8 

59.4 

81.8 

Post-SNe 

328.5 

319.2 

107.2 

126.6 

5 Myr 

64.6 

57.6 

54.1 

64.3 

25 Myr 

84.9* 

67.8* 

61.8* 

62.4 


Bold: < 50% retained 

b 0 Myr = tsNe? =6.1 Myr from start of simulation 
Pre-SNe 50kyr prior to SNe, Post-SNe, lOkyr after SNe 
*: Increased Mass after some infall 
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radii, Vs and Tvir, are included for each model. The M60, 
M65 and M70 models retain essentially all of the gas 
within Tvir for the cooling models, while the M55 case 
(not shown) is evacuated completely. For the adiabatic 
cases (not shown), the M60 halo is evacuated completely 
and the M65 halo appears to lose more than half the 
gas after 30 Myr. In this instance, as the gas reacts 
to the radiation field, it becomes photoheated, expands, 
smooths out and drops in density as a result. The num¬ 
ber of recombinations locally declines and the radiation 
propagates into more of the surrounding gas. 

In Figs. 12 to 17, we observe a number of important 
facts about the pre-ionization phase. For the centred 
models, the Stromgren radius of the HII region falls out¬ 
side (M60), onto (M65) and inside (M70) of the scale 
radius. From the zoomed figures, we observe how the 
photoheating smooths out the substructure. The effect 
is so strong that we observe a cellular structure as the ion¬ 
ization fronts of many expanding knots overlap, which is 
particularly clear in the colour composite figure (Fig. 18). 
In the same figure, we see broad channels appear where 
the radiation is able to escape to the outer halo over a 
substantial solid angle. Even before the SN event, the 
expanding Stromgren sphere has produced a diffuse neu¬ 
tral shell of gas around the hot bubble cavity. In the off- 
centred models, the expanding Stromgren sphere com¬ 
presses the core region gas and the radiation is blocked 
from propagating over at least 27r sr. This becomes im¬ 
portant in the post-supernova phase because substantial 
gas begins to re-accrete to the centre along this axis. 

5.2. Supernova phase 

Our complete set of simulations (see Table 3) comprises 
24 cases for supernova events without pre-ionization with 
the massive progenitor, and 6 cases where the ionization 
of the supernova progenitor is included. All 30 cases 
were run at LOO and LOl, with a subset of the cooling 
models run at L02 to allow adequate sampling of the 
developing supernova shell (see Fig. 5). L02 is always 
found to be adequate to resolve the metal-enriched shell 
which starts as a diffuse hot envelope which becomes 
thinner under expansion. The models do not suffer the 
overcooling through inadequate resolution identified by 
Dalla Vecchia & Schaye ( 2012 ). 

We discuss in detail only the clumpy, cooling models 
for the centred and off-centred supernova cases where 
the progenitor ionization is included. We present these 
results for three different halo masses (MOO, M65, M70)^. 
Tables 5 and 6 include the total gas mass and energy 
retained within r = r^ and r = Tvir after periods of t = 
5 Myr and t = 25 Myr. The halos that retain less than 
50% of their original mass or energy have their retention 
fractions emboldened. 

Each of the models (MOO, M05, M70) are discussed in 
more detail below. In Eigs. 12 to 17, we present time- 
dependent 2D slices of the gas distribution before and 
after the SN explosion. These are not summed along the 
third dimension and so therefore do not illustrate the full 
extent of the ionized gas (cf. Eig. 18). Unless otherwise 
stated, all computations are performed over a 048^ grid 
with two levels of the hierarchy (LOO, LOl) illustrated in 

^ All simulations can be viewed at 
http://miocene.anu.edu.au/smallgalaxy . 


Eig. 5 (see §3.4). In a subset of cases, it was necessary 
to investigate a 3-level hierarchy (LOO, LOl, L 02 ) over a 
19443 which we clarify below. 

It is worthwhile to pay particular attention to Eig. 19. 
Here we present ID plots of the evolution of the gas mass 
(A-E) and energy (G-L) retention fraction over much 
longer timescales than shown in the half-tone figures (up 
to 25 Myr). In each plot, we present two sets of curves 
for retention within r <rvir (upper) and r <r^ (lower). 

The left panels (A, G, E, G, I, K) show all combinations 
— clumpy (solid curves) vs. smooth (dashed curves), off- 
centred (light curves) vs. centred (heavy curves), adia¬ 
batic (red curves) vs. cooling (black curves). These plots 
are particularly illustrative and allow us to discuss the 
role of internal and kinetic energy in the developing ex¬ 
plosion. Eor all models, the adiabatic cases (red curves) 
confirm earlier work in that there is little or no gas reten¬ 
tion (Mac Low & Eerrara 1999). In the right panels (B, 
D, E, H, J, L), we show only the cooling clumpy models 
for the centred and off-centred cases; the gas retention is 
clearly more favourable for the larger halo masses. 

One important insight emerges which explains why we 
have arranged the figures in this way. In each row, we 
compare two distinct halo masses (e.g. M55 vs. M60 
in the top row) to illustrate the impact of pre-ionization 
due to the supernova progenitor. The left panels (e.g. 
A) have no pre-ionization in contrast to the right pan¬ 
els (e.g. B) which include this effect. Broadly speaking, 
the results on retention are similar. The same holds true 
when we compare M60 models without pre-ionization to 
M65 models with pre-ionization, and M65 models with¬ 
out pre-ionization to M70 models with pre-ionization. 
These are important limiting cases when one considers 
enrichment in halos that do not have an important pre¬ 
ionization phase, e.g. for star formation histories limited 
to low-mass and/or Type la supernovae. 

5.3. M60 models 

The time-dependent behaviours of the M60 central and 
off-centred explosions are illustrated in Eigs. 12 & 13 re¬ 
spectively. The evolutionary phases of an explosive wind 
are well understood analytically (Weaver et al. 1977) al¬ 
though we see additional substructure due to the com¬ 
plexity of the initial medium. The supernova ignites 
within a fully ionised, diffuse nebula confined by a neu¬ 
tral gas shell. A dense radiative shell quickly forms and 
drives a shock into the neutral shell. The momentum- 
driven shell does not stall and essentially evacuates all 
gas within the halo after 15 Myr. In the off-centred case, 
one hemisphere is initially evacuated within 15 Myr, with 
the opposite hemisphere losing most but not all of its gas. 
At late times, there is some re-accretion from the far side 
leading to < 10 % gas retention within the core regions. 
The energy input from the pre-ionization phase, in com¬ 
bination with the SN explosive energy, conspire to doom 
the M60 halos to invisibility. 

Numerieal artefaets. In order to examine the impact of 
artefacts, this case is treated with a 3-level grid. No 
such anomalies are seen in the pre-ionization phase. In 
the SN event, numerical artefacts evident in the first 600 
kyr near the origin (Eig. 20) are caused by the carbuncle 
instability (cf. Sutherland 2010). This effect is ampli¬ 
fied by strong cooling in the early stages (< 125 kyr) 
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Fig. 12 .— M60 photoionization (centred source) during the pre-supernova and supernova phases at four different epochs (0.5, 2.0, 6.0, 
6.5 Myr). The left and right figures in each panel show the projected neutral and ionised gas column respectively. The supernova ignites at 
6.0 Myr. The detailed geometry seen here shows how radiation breaks out from the inner core by 0.5 Myr along certain directions although 
the timescale is longer (2.5 Myr) overall. The density variations make the ionization fronts complex due to shadowing by the densest gas. 
However the dense gas is too sparse to present a significant barrier. Immediately prior to the SN explosion, the HII region is fully ionised, 
heated and highly smoothed the gas, which improves the SNe - ISM coupling, further increasing the ability of the SNe to clear the halo. 
The warm gas is less bound than initially as well, and by 6.5 Myr the SNR has cleared the core region completely. 
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Fig. 13.— M60 photoionization (off-centred source) during the pre-supernova and supernova phases at four different epochs (see Fig. 12). 
The supernova ignites at 6.0 Myr. The escape of photons is highly asymmetric and ioization is complete (except on the diametrically 
opposite side of the core) by 2 Myr. By 6 Myr, the core is heated, ionised, and smoothed, so the SN is quite efficient in clearing the halo. 
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Fig. 14.— M65 photoionization (centred source) during the pre-supernova and supernova phases at four different epochs (see Fig. 12). 
The supernova ignites at 6.0 Myr. The non-uniform gas presents low density channels, so that by 2.0 Myr, the spiky ionization fronts have 
escaped the core in the central model. By 6.0 Myr, the entire core region is photoionised. 
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Fig. 15.— M65 photoionization (off-centred source) during the pre-supernova phase at four different epochs (see Fig. 12). The supernova 
ignites at 6.0 Myr. The off-centred Hll region runs out down the halo density gradient rapidly on the side near the star, but struggles to 
ionise through the core to the opposite side. The geometry allows the ionization to affect a larger volume in the outer halo than indicated 
by the integrated analysis in Fig. 7. The SN occurs again in a smoothed region and reaches a significant fraction of the core radius by 6.5 
Myr. 
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Fig. 16.— M70 photoionization (centred source) during the pre-supernova and supernova phases at four different epochs (see Fig. 12). 
The supernova ignites at 6.0 Myr. The photon ffux is insufficient to ionise a significant part of the core volume, and the Hll region is stiffed. 
The SN occurs in a similar environment to the non-preionised models and the outcomes are similar. Essentially all of the gas survives the 
impact of the SN event. 
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Fig. 17.— M70 photoionization (off-centred source) during the pre-supernova and supernova phases at four different epochs (see Fig. 12). 
The supernova ignites at 6.0 Myr. The off-centred Hll region is able to reach the edge of the core, and the falling density slope to the outer 
halo ensures that the photoionization fronts accelerate outward on the star side. However a large volume of the halo remains unaffected. 
The SN is able to just reach the core radius by 6.5Myr, but has only a small effect on the bulk off the halo. The SN occurs in a similar 
environment to the non-preionised models and the outcomes are similar. Almost all of the gas survives the impact of the SN event. 
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Fig. 18.— The development of the pre-ionization and SN phases for the M60, M65 and M70 cases for four different epochs; the top three 
figures are for a central source and the lower three figures are the corresponding off-centred cases. The snapshots shows the distribution of 
ionised and neutral gas before, during and after the SN event at 6.0 Myr. This figure is a colour composite of the information contained 
in Figs. 12 to 17. Each image is a density slice through the centre of the halo. The red shows the distribution of fully ionised gas; the blue 
is the neutral gas; mauve corresponds to partially ionised gas. 
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and leads to the ‘clover leaf’ pattern at 6.3 Myr inside 
Vs (LOl). Our algorithm stabilises the instability when 
the radiative thin shell becomes only partially resolved 
(< 3 cells thick) at later times, such that any structure 
that evolves in the interaction between the spherical shell 
and the cartesian grid remains at all later stages. The 
instability is strongest where the radiative shock front 
becomes orthogonal to the underlying grid, in this case 
along the cartesian axes, giving rise to the boxy struc¬ 
tures with apices along the x— and y— axes seen here 
at 16 Myr (LOO). After the local density has dropped so 
that the blast-wave is no longer radiative, the structure 
then simply expands with the now adiabatic blast. The 
additional L2 level is needed to better resolve the early 
blast, and the remaining visible non-spherical structure 
has radius errors of approximately 10 — 15%. 

Gas retention and energeties. We now discuss the issue 
of gas retention in the MOO models in the context of the 
energy (internal+kinetic) of the disturbed gas. Table 6 
shows the amount of energy retained at critical radii, 
r = Vs and r = Tvir, immediately before and after the 
supernova explodes, and at 5 Myr and 25 Myr after the 
supernova event. We show density slices through our 
clumpy MOO models over this time frame in Figs. 20 and 
21 . 

We recall that during the pre-ionization phase, radia¬ 
tion leaks out from the centre in the first few Myr leading 
to a fully developed, low density nebula by 6 Myr (see 
§3.5). In Fig. 20, the central SN explodes into a warm 
low-density nebula, superheats the gas through shock 
heating, and removes essentially all of the gas from the 
dark matter halo out to Tvir after 25 Myr. This is borne 
out in Table 5. 

For the off-centred case in Fig. 21, the progenitor su¬ 
pernova is placed at the half gas-mass radius to reflect 
the fact that most stars form off-centre. During the pre¬ 
ionization phase, radiation leaks out from the centre in 
the first few Myr (see §3.5), and a conic ‘zone of avoid¬ 
ance’ is set up diametrically opposite the progenitor. By 
6 Myr, a low density nebula extends across most of the 
core region. The SN explodes into an off-centred, a warm 
low-density nebula, superheats the gas through shock 
heating, and removes much of the gas asymmetrically 
from the dark matter halo out to Tvir after 25 Myr. The 
early zone of avoidance leads to far-side gas re-accreting 
into the core region after the SN event has passed, but 
less than 10% of the original gas remains within Tvir on 
long timescales. 

5.4. M65 models 

We now discuss gas retention in the M65 models within 
the context of the mass and energy of the disturbed gas. 
The time-dependent behaviours of the M65 central and 
off-centred explosions are illustrated in Figs. 14 & 15 
respectively before and ineluding the explosion, and in 
Figs. 22 & 23 respectively after the supernova event. To 
aid comparison, the projected neutral and ionised gas in 
Figs. 14 & 15 are combined into a colour composite in 
Fig. 18. 

Tables 5 and 6 show the percentage of mass and energy 
retained at and Tvir just before and after the supernova 
explosion, and after delays of 5 and 25 Myr. We see that 
more gas mass is retained in the core in the off-centred 


case compared to the central case after 25 Myr post su¬ 
pernova. The tabulated energy retentions are percent¬ 
ages relative to the binding energy of the gas bound by 
the dark halo. A system loses energy either through cool¬ 
ing (e.g. radiative, adiabatic) or mass loss; conversely, a 
system gains energy either through heating (e.g. shocks, 
UV radiation) or accretion. In the off-centred model, we 
see the effects of mass loss through the virial radius at 
the same time that the core is accreting gas. 

In Fig. 14, the supernova ignites at 6.1 Myr within 
a fully ionised, diffuse nebula confined by a neutral gas 
shell. A dense radiative shell quickly forms and drives a 
shock into the neutral shell. A lot of energy is radiated 
away in the time frame 0.2 to 2.5 Myr. The momentum- 
driven shell stalls at 1 Myr. At this point, the shock slows 
down and becomes non-radiative. The drop in temper¬ 
ature leads to a rapid increase in density. The neutral 
shell then becomes energy-driven and more diffuse as it 
propagates outwards at about the halo escape speed and 
ultimately stalls at the virial radius. Most of the gas just 
survives in the M65 progenitor ionization + SN model. 

In Fig. 15, the off-centred HII region runs out down the 
halo density gradient rapidly on the side near the star, 
but struggles to ionise through the core to the opposite 
side. The geometry allows the ionization to affect a larger 
volume in the outer halo than indicated by the integrated 
analysis in Fig. 7. The SN occurs again in a smoothed 
region and reaches a significant fraction of the core radius 
by 6.5 Myr. 

In these models, no more than half the gas survives 
the SN explosion, is enriched with the specific abundance 
yields of the discrete SN event. Our highest resolution 
simulations reveal why cooling is so effective in retaining 
gas compared to other factors. In the early stages, the 
super-hot metal-enriched SN ejecta exhibit strong cool¬ 
ing, leading to much of the explosive energy being lost. 
This remains true regardless of the assumed metallicity 
floor at the time of the explosion. 

In Figs. 22 & 23, we explore the post-supernova phase 
over larger radial scales (LOO) and on longer timescales. 
The UV radiation from the supernova progenitor escapes 
along some sight lines over the full spatial extent of the 
models. In both cases, the development of the supernova- 
driven bubble is well defined. A thick HI shell forms 
by 10 Myr ahead of the thin metal-enriched shell and 
this material starts to reaccrete to the core by 25 Myr. 
In the off-centred case, the central column densities are 
higher because the shock front wraps around the core 
and compresses it. 

5.5. M70 models 

We now consider the retention of gas and energy in the 
M70 models. The HI and HII column densities before 
and just after the supernova are shown in Figs. 16 and 
17, while Figs. 24 and 25 show the HI density evolution 
after the supernova. For ease of comparsion. Figs. 16 
and 17 are also shown on the colour composite figure 
Fig. 18. The M65 models walk the tightrope between 
survival and disruption, retaining half their baryons at 
best. Tables 5 and 6 show percentages of gas and energy 
retained just before and after the explosion, as well as 
after 5 and 25 Myr. The M70 models clearly survive 
a single supernova event, losing very little gas and by 
25 Myr contain more gas than in their initial state as a 



Lowest Mass Galaxies 


23 


result of accretion. As a result of cooling, the energies 
within both and Tvir are much lower than the initial 
values at both 5 Myr and 25 Myr, suggesting that the 
gas is not flowing outwards. 

For the M70 models, for the higher total gas mass, 
more than a single event could explode at a given epoch. 
But this is a complication that we do not consider here 
because the combined impact of the two events will de¬ 
pend on their masses (one is likely to be in the low-mass 
limit) and locations within the halo. In Fig. 16, the cen¬ 
tral supernova ignites at 6.0 Myr. The photon flux is 
insufficient to ionise a significant part of the core vol¬ 
ume, and the HII region is stifled. In the off-centred case 
in Fig. 17, the off-centred HII region is able to reach the 
edge of the core, and the falling density slope to the outer 
halo ensures that the photoionization fronts accelerate 
outward on the side containing the star. The supernova 
is able to just reach the scale radius by 6.5 Myr, but has 
only a small effect on the bulk of the halo. 

The post-supernova phase is shown in Figs. 24 and 
25. For the central supernova model, the bubble at 
t = 2.5 Myr is smaller than in the M65 case and sweeps 
up denser gas. By 10 Myr the bubble collapses and dense 
gas accretes back to the centre. In the off-centred case, 
the behaviour is similar, although a small amount of ra¬ 
diation reaches the edge of the grid on the side nearest to 
the supernova, while the gas on the opposite side of the 
halo is largely unaffected. As in the M65 off-centred case, 
the central densities are increased at 2.5 Myr when the 
shock compresses the core, however the effect is weaker 
than in the M65 model because the shock is unable to 
wrap around the core. 

Our models consider dwarf galaxies to be isolated, but 
in practice the vast majority reside in the extended halo 
of a host galaxy, albeit beyond the virial radius. In the 
M65 and M70 models, there is clear evidence of the ex¬ 
pelled gas beginning to fall back after 20 Myr. A proper 
consideration of this gas requires that we consider the 
environment of the halo which is beyond the scope of the 
current work. We leave this to a later study. 

5.6. Summary 

As described in the preceding sections, the summary 
plots for the single supernova event models in Fig. 19 
(supported by Figs. 12 to 17) reveal that: 

1. For all M55 models and M60 adiabatic models, es¬ 
sentially none of their original gas mass remains 
within a scale radius at 25 Myr. Thus subsequent 
star formation does not occur. 

2. For all M60 cooling models, only a few percent of 
the original gas mass remains within a scale radius 
at 25 Myr. No significant star formation can occur. 

3. For the M65 models, all of the cooling scenarios 
retain 70-80% of the original gas inside the virial 
radius at 25 Myr. This is sufficient for ongoing 
star formation. None of the adiabatic models retain 
more than 20% of the original gas. 

4. In the borderline M65 cooling cases, off-centred 
models retain 2-3 times more gas within a scale ra¬ 
dius but the supernova metals mix less efficiently 
with this gas. 


5. For the M70 models, all of the gas survives a su¬ 
pernova event and its radiative progenitor. We 
stress that we have not considered clustered multi¬ 
supernova events which can occur in these models. 

6. Cooling is the most important factor governing gas 
retention within the virial radius at late times. 

7. Clumpiness is far more important in the adiabatic 
case than in the cooling case for gas retention. 

6. DISCUSSION 
6.1. Background 

The least massive galaxies lie at the bottom of the 
CDM hierarchy. Our focus here is not their role as build¬ 
ing blocks in galaxy evolution, but what they can tell 
us about the yields of the first stars prior to the reion¬ 
ization epoch, and the following generations during and 
immediately after reionization (e.g. Norris et al. 2010; 
Bland-Hawthorn et al. 2010; Frebel & Bromm 2012; Rit¬ 
ter et al. 2014). At least some of the first galaxies must 
have resided within low mass halos. We anticipate that 
such systems, if they can be identified in the present day 
universe, are likely to carry the most ancient chemical 
signatures because of their age and because of the rela¬ 
tively low number of enrichment events over cosmic time. 

If we can identify such objects, it should be possible 
to obtain a relatively clean signature of the primordial 
yields from the earliest stellar generations. The star- 
formation efficiency is greatly diminished in low-mass 
galaxies, thereby reducing the total number of enrich¬ 
ment events. There are likely to be instances where a 
large fraction of the metals are blown out of the galaxy, 
and this may obscure the signatures of earlier genera¬ 
tions of stars. This can occur when the relatively rare 
supernovae that go off in these galaxies are correlated in 
time. But we would still expect to find the signatures of 
a single enrichment event in at least some intrinsically 
low-mass galaxies. 

We now discuss our models in light of two recent dis¬ 
coveries: ultra-faint dwarf galaxies (near field) identified 
in all-sky imaging surveys (Simon & Geha 2007; Kirby 
et al. 2008), and very metal poor damped Lya systems 
(far field) identified along quasar sight lines (Erni et al. 
2006; Cooke et al. 2010). The two classes of objects in 
turn although a more extensive discussion is given else¬ 
where (Webster et al. 2014, 2015a,b). We refer the reader 
to Figs. 1 and 2 that show slow evolution is expected be¬ 
tween the formation redshift at z = 10, the epoch of the 
DBAs {z = 2 — 4), and present-day ultra-faint dwarfs. 
The overall effect is additional dark matter that is mostly 
accreted to the outer envelope of the dwarf, rather than 
onto the central cusp (Fig. 4). 

6.2. Ultra-faint dwarfs 

In recent years, the Sloan Digital Sky Survey has dou¬ 
bled the number of known dwarf galaxies that orbit the 
Galaxy (for a review, see Tolstoy et al. (2009)), particu¬ 
larly at the low mass end. These newly discovered dark- 
matter dominated galaxies are uniformly very old and 
cover a wide range in luminosity down to ^ 10^ Lq. Un¬ 
like globular clusters, they show a wide range of eccen¬ 
tricities which may indicate that some are undergoing 
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tidal disruption (Martin et al. 2008). Objects less lu¬ 
minous than 10^ Lq are referred to as ultra-faint dwarfs 
(UFDs). 

There has been extensive discussion concerning the in¬ 
tegrated masses of UFDs as observed today (Simon & 
Geha 2007; Strigari et al. 2008; Martin et al. 2008; Bovill 
& Ricotti 2009; Wolf et al. 2010; Frebel et al. 2014). Some 
studies quote a characteristic UFD mass within 300 pc, 
i.e. M 300 ~ 10^ Mq (Strigari et al. 2008; Walker et al. 
2009) although this is an arbitrary physical scale that is 
not easily related to well known scaling relations. The 
characteristic virial mass Mvir today is highly uncertain 
with estimates ranging from 10^ to 10^ Mq. For a dark 
halo described by the NFW profile, the mass within a 
fixed radius gives the mass within any other radius, as¬ 
suming this is an accurate description of the underlying 
gravitational potential. For our work, we need to know 
the virial mass when the UFD first formed. The halo 
parameters in cosmological simulations have a strong de¬ 
pendence with redshift making any estimates rather un¬ 
certain (Power et al. 2003; Duffy et al. 2008). We have 
adopted Mvir ^ 10^ Mq because the anticipated total 
population of UFDs may be too large to support larger 
progenitor masses. 

The simulated ‘observables’ from our models are 
broadly consistent with the most metal-poor UFDs 
(Webster et al. 2014, 2015a). The data for these systems 
are summarised by Frebel & Bromm (2012, Fig. 1). The 
metallicity range for our models is -1.0 < [Fe/H] < -4.0 
with a mean [a/Fe] enhanced by a factor of a few with 
respect to the Solar value depending on the initial mass 
function at the time of formation (Woosley & Weaver 
1995). We have begun to consider the r and 5 process 
signatures in the gas, work that is to be presented in later 
papers. The internal kinematics of the baryons are less 
than 4 km s“^ within the inner scale radius, consistent 
with observations (Walker et al. 2009). 

Brown et al. (2014) find that the star formation history 
of UFDs is marked by one or two bursts of star formation 
before truncation by a global event. By construction, 
such a history is consistent with our models. UFDs have 
stellar populations with mean metallicities in the range 
-2.5 < [Fe/H] < -2.0 (Kirby et al. 2013). Within the 
framework of our models, this requires (i) 1 — 2 SN events 
in a low mass halo, specifically the M65 event-centred 
models (§5.4); (ii) 3 — 10 SN events in a starburst within 
a dark matter halo an order of magnitude more massive 
(e.g. M70). These cases can give rise to measurably 

distinct chemical histories, star formation histories and 
star cluster populations. For example, [a/Fe] is enhanced 
in starbursts relative to quiescent star forming regions. 
Direct evidence for this comes from x-ray observations of 
starburst winds (Martin et al. 2002). 

At the time of publication, we note that Ohorbe et al. 
(2015) have presented a new model for UFDs where the 
baryons reside within dark matter halos that are orders 
of magnitude more massive than we consider (Mvir~ 
3 X 10^ Mq). Given that all of the baryons are retained 
within our M70 models, for the same initial gas fraction, 
we cannot easily distinguish their model from ours in 
terms of the predicted chemical abundances. 

In summary, if UFD masses are largely unchanged 
since their early formation, we would expect these halos 


to show significant scatter in [Fe/H] at low metallicity, 
and for some of the oldest stars to constrain the yields 
of early stellar generations (e.g. Frebel & Bromm 2012). 

6.3. Very metal poor DLAs 

An interesting development is the recent discovery of 
very metal-poor damped Lya systems (DLA) along QSO 
sight lines (Pettini et al. 2008; Penprase et al. 2010). 
If these are protogalactic structures that have recently 
formed from the IGM, there is the prospect of identifying 
the chemical imprint of early generations of stars (e.g. 
Pettini et al. 2002). Interestingly, two of these systems 
with [Fe/H] ^ -3 may bear the hallmarks of early stellar 
enrichment (Erni et al. 2006; Gooke et al. 2010, 2015). 

First, Erni et al. (2006) identify a DLA towards the 
QSO Q0913+072 (z = 2.785) with an iron abundance 
characteristic of the IGM at that redshift. The G, N, O, 
Al, Si abundances show an odd-even effect reminiscent 
of the most metal-poor stars in the Galactic halo. This 
pattern is created in models where the neutron flux is 
low (e.g. Heger & Woosley 2002), presumably due to the 
low overall metal abundance. A more striking signature 
is the strong [N/H] depletion which Pettini et al. (2002) 
has argued is further evidence for a system that has re¬ 
cently formed from the IGM. Erni et al. (2006) argue 
that the abundances appear to agree with 10—50 Mq 
zero metallicity Pop HI models. 

Secondly, Gooke et al. (2015) identify half a dozen DLA 
objects that resemble the properties of UEDs with mean 
[Ee/H] in the range -2.5 < [Ee/H] < -2.0 and line-of-sight 
dispersion velocities of 1 < (Jlos < 5 km s“^ corrected 
for thermal and turbulent broadening. These physical 
parameters are in line with our M65 and M70 models 
(Webster et al. 2015b). 

Einally, we examined whether the line-of-sight proper¬ 
ties of our models would quality them as DLAs. Eor all 
very metal poor DLAS, Gooke et al. (2015) infer the total 
mass of neutral gas to be 10^ < Mgas ^ 10^ Mq within 
a linear scale of < 30 — 1300 pc although their UED can¬ 
didates are at the low end of these ranges. There are 
no sight lines in the M55 models, and few sight lines in 
the M60 models, that reach the DLA threshold. But the 
M65 and M70 models fare better with a significant frac¬ 
tion of sight lines within a radius of ^ being classified 
as a DLA-type spectrum. Eurthermore, the inferred gas 
metallicity of [Ee/H] ^ -3 falls well within our predicted 
range. This gas is mostly self-shielded from any exter¬ 
nal radiation and therefore does not require an ionization 
correction. 

7. GONGLUSIONS 

To investigate the issue of gas retention in the lowest 
mass dark-matter halos, we performed a series of exper¬ 
iments in the form of high-resolution radiative hydrody¬ 
namic simulations dominated by the gravitational poten¬ 
tial of the dark matter. We find that the broad properties 
of UEDs and very metal poor DLAs are consistent with 
the host halos having intrinsically low mass. Our calcula¬ 
tions are based on the critical period a few Myr before the 
SN explodes, and up to 25 Myr after the event. These 
cases are explored over much longer timescales (^ 600 
Myr) elsewhere which allows us to consider the effect of 
subsequent SN events and the homogenising effects of re- 
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accreting gas dislodged by earlier events (Webster et al. 
2014, 2015a,b). 

In summary, dark matter halos with virial masses as 
little as 3 X 10^ Mq are able to retain some gas, a limit 
that is an order of magnitude lower than found previ¬ 
ously (Mac Low & Ferrara 1999). The effects of cooling, 
the dumpiness of the medium, and the position of the ex¬ 
plosion were studied. Consistent with earlier work, the 
baryons do not survive in smooth adiabatic or smooth 
cooling models in the event of a supernova. Models with 
radiative cooling, clumpy media and off-centred explo¬ 
sions are most favourable to gas retention. Our highest 
resolution simulations reveal why cooling is so effective in 
retaining gas compared to any other factor. In the early 
stages, the super-hot metal-enriched SN ejecta exhibit 
strong cooling, leading to much of the explosive energy 
being lost. This remains true regardless of the assumed 
metallicity floor at the time of the explosion. 

An added complication is to consider asymmetric ex¬ 
plosions in core-collapse supernovae (Grefenstette et al. 
2014) rather than the isotropic explosions considered 
here. But we believe these do not lead to radically differ¬ 
ent outcomes when compared to the off-centred, isotropic 
supernova models. To first order, the overall mixing is 
dictated by geometry such that an energetic pulse experi¬ 
encing any asymmetry leads to broadly similar metallic¬ 
ity spreads. If confirmed by later work, this principle will 
greatly reduce the parameter space of future simulations. 
It remains to be seen if asymmetric vs. symmetric explo¬ 
sions lead to the ejection of different abundance yields. 

We recognise at the low mass limit explored in our 
models, the predictions are subject to stochastic effects 
due to the rarity of supernova events. Any assumed ini¬ 
tial mass function will be poorly sampled in this limit. In 
rare instances, simultaneous supernova events are statis¬ 
tically possible but these are more likely in higher mass 
halos than considered here. Gonversely, it is conceivable 
that low-mass star formation can proceed at the adopted 
metallicity threshold of [Fe/H]=-4 (Bromm & Yoshida 
2011) such that no prior supernova activity is required 


beyond the enrichment from the first stars. This leads 
to a small fraction of stars extending down to the metal¬ 
licity floor, a signature observed in only a few cases, e.g. 
Bootes I (Gilmore et al. 2013) and Segue I (Norris et al. 
2010 ). 

This work has identified new lines of enquiry. In par¬ 
ticular, the central gas densities are chosen to limit the 
gas mass inside one scale radius to about 10% of the 
dark matter mass on the same scale (§2). Such a highly 
concentrated baryon fraction is to be expected because 
the star formation efficiency is likely to be very low be¬ 
fore the onset of metal (and therefore dust) production 
(Krumholz & Dekel 2012). In future papers, we will con¬ 
sider higher baryon fractions including those that dom¬ 
inate over the dark matter core. This requires us to in¬ 
clude self gravity into the Fyris code through the use 
of fast fourier transforms over the refined adaptive grid. 
In these models, we treat the pre-supernova phase with 
higher resolution grids. We have already begun to imple¬ 
ment self-consistent cooling arising from variable enrich¬ 
ment over the gas distribution. Our future models will 
include the formation of dust and molecules. A complete 
understanding of UFDs will need to include the details of 
the subhalo orbit within the evolving Galactic potential 
and the coronal gas distribution responsible for the gas 
stripping (Nichols et al. 2012, 2015). 
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Fig. 19.— Gas mass (A-F) and energy (G-L) retention inside the M55, M60, M65 and M70 models. Panels A, G, E and G, I, K are 
for SN explosions without pre-ionization (Type la SN); panels B, D, F and H, J, L are for SN explosions with pre-ionization (Type II 
SN). For the former, we start the models at the explosion so there is no data prior to the event. The time of the SN event is shown as a 
vertical dashed line; the impact of pre-ionization (see Table 5) is clearly evident. The key to all curves is given in both figures. The top 
set of curves in each panel is for quantities within a radius the bottom set is for quantities within rg. The mass fraction is defined 
relative to the original gas mass within the scale radii as indicated at t = 0. We show both the adiabatic (red curves) and cooling models 
(black curves) for both clumpy and smooth media. Gonsistent with earlier work, the smooth adiabatic case is much lossier than the case 
for cooling halos. Note that the adiabatic case is energy conserving until gas is lost over a boundary. In the energy retention figures, note 
that M70’s explosive energy has only a marginal impact on the energy retained from the progenitor phase. 
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APPENDIX 

APPENDIX A - ABUNDANCES 

The host halo gas metallicity was taken to be [Ee/H] = -4.0. The abundances, by number relative to hydrogen 
= 1.0, are shown in Table 7. The solar abundances are based on Asplund et al. (2005), with a mild alpha element 
enhancement pattern: Carbon (xl.5), Oxygen(x5.0), Magnesmm(x2.5), Sihcon(x2.5), and Calcmm(x2.5) (K13). 
Helium evolution follows (Dopita et al. 2006), uhcI'^h = 0.0737+0.0024 Z/Zq. 

The supernova ejecta composition was modelled on the well studied M25 supernova, cf. Kobayashi et al. (2006); 
Nomoto et al. (2006). Using the Kobayashi compilation yields 2.79 My) of oxygen in the ^ 4.4 My) of total metals 
ejected, with 0.072Me of Fe 

Table 7 shows the SNe mass yields, and the abundances by number, for the ejecta and ambient medium, with a final 
column showing solar abundances ratios for reference. 

APPENDIX B - ERACTAL INTERSTELLAR MEDIUM 
Log-normal density distribution 

We use a log-normal distribution to describe the single-point statistics of the density field of our non-uniform ISM. 
The log-normal distribution is a skewed continuous probability distribution. Unlike the normal distribution, it has a 
non-zero skewness, variable kurtosis, and in general the mode, median and mean are unequal. 

The log-normal distribution appears to be a nearly universal property of isothermal turbulent media in experimental, 
numerical and analytical studies (e.g. Nordlund & Padoan (1999); Warhaft (2000); Pumir (1994)). Moreover, it is 
encouraging that the log-normal distribution is the limiting distribution for the product of random increments, in the 
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TABLE 7 

Gas Phase Abundances, 

[Fe/H] = -4.0, Neutral: /j, = 1.21, Ionised: ii = 0.59 



SN25 Z=0 

SN25 Z=0 

‘Primordial’ 

Solar 

Element 

Mass Yields 

Abundances 

Abundances 

Ratios 

H(X) 

10.60 

0.00 

0.00 

0.00 

He (Y) 

8.03 

-0.72 

-1.13* 

-1.01 

Metals(Z) 

4.41 

-1.62 

-6.27 

-3.09 

C 

2.94E-01 

-2.63 

-7.28 

-4.10 

N 

5.91E-04 

-5.40 

-10.05 

-4.79 

O 

2.79E+00 

-1.78 

-6.43 

-3.34 

Ne 

5.33E-01 

-2.60 

-7.25 

-3.91 

Na 

1.03E-03 

-5.37 

-10.02 

-5.75 

Mg 

1.20E-01 

-3.33 

-7.98 

-4.42 

A1 

8.08E-04 

-5.55 

-10.19 

-5.61 

Si 

3.51E-01 

-2.92 

-7.57 

-4.49 

S 

1.86E-01 

-3.26 

-7.91 

-4.79 

Ar 

3.14E-02 

-4.13 

-8.78 

-5.20 

Ca 

2.48E-02 

-4.23 

-8.88 

-5.64 

Fe 

7.38E-02 

-3.90 

-8.55 

-4.55 

Ni 

4.42E-04 

-6.15 

-10.79 

-5.68 


* Helium from BB NS, metals scaled to [Fe/H] = -4.0. 


same way that the normal distribution plays that role for additive random increments. It is thus compatible, at least 
conceptually, with a generic cascading process consisting of repeated folding and stretching. 

With a log-normal distribution, which is on average isotropic, the natural logarithm of the ISM density field is a 
Gaussian which has a mean m and variance . The probability density function for the log-normal distribution of the 
mass density p is. 


p{p) = 



— (In p — m)^ 
2^2 


(Bl) 


The mean p and variance of the density are given by 


/i = exp[m + 5^/2]. (B2) 

ap=p‘^ (exp[5^] - 1). (B3) 


In these simulations we adopt // = 1.0, = 5.0, as our standard log-normal distribution. With these parameters, 

densities below the mean comprise one quarter of the mass, and occupy three quarters of the volume, and the mean 
is approximately 20 times the mode. See Sutherland & Bicknell (2007) for more details. 


Power-law density strueture 

The two-point structure of a homogeneous turbulent medium may be described in Fourier space. We denote the 
Fourier transform of the density p(r) by F(k) (where k is the wavenumber vector). The isotropic power spectrum 
D{k) is the integral over solid angle in Fourier space of the spectral density F(k)F*(k). In three dimensions: 

D{k) = J k^F{k)F*{k)dn. (B4) 

. Even if the spectral density is anisotropic, the angular integral averages the spectral density into a one dimensional 
function of k only. For a power-law dependence on k, D{k.) oc k“^ and p = 5/3, the spectrum is referred to as 
Kolmogorov turbulence. 

In order to simultaneously achieve log-normal single-point statistics and a power-law self-similar structure, we have 
implemented the practical method developed by atmospheric scientists, for constructing two and three dimensional 
terrestrial cloud models, which are used in radiative transfer calculations (Lewis & Austin 2002). 

The remaining choice in this procedure is to select the range of wave numbers over which to generate the fractal, in 
particular the minimum wave number kmin, which determines the largest structure scale in the resulting fractal with 
respect to the spatial grid. The main domain of interest on the intermediate level in the simulations covers at least 
two dark matter scale radii in each model. We choose kmin = 8 over this domain, corresponding to approximately two 
cloud structures per scale radius. On the upper level, the domain covers three times the extent, the fractal is extended 
in a kmin = 24 model. The lowest level is initially populated with a simple resampled version of the intermediate level 
density field. Each level of the simulations is sampled in a 216 ^ cells and so even the coarse kmin = 24 upper level has 
well resolved small clouds. 
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Fig. 20.— M60 Central Clumpy explosion commencing at the end of the progenitor ionization (pre-ionization) phase. The SN explodes 
into a warm low-density nebula pre-heated by the progenitor. The SN wind superheats the gas through shock heating, and removes 
essentially all of the gas from the dark matter halo out to rvir after 30 Myr. For the uppermost 4 panels (LOO), the inner circle has radius 
rg; the outer circle has radius rvir- Each panel has been magnified in the next 4 panels (LOl); the circle has radius rg. All images are 
central plane density slices at z = 0.0. Log of the normalised density variable in units of 1.0 x 10“^^ g cm“^: greyscale minimum = —3.0 
(white), maximum = 2.0 (black). A scale bar in parsecs is shown. 


Equilibrium turbulent distribution 

Here we are considering purely spherical distributions, neglecting rotational formation of disk structures. The reasons 
are two-fold. First, flattening due to rotation will provide a natural axis for the escape of the SN energy and potentially 
greatly reduce the efficiency of the supernova in evacuating the region of potentially star-forming gas. We are looking 
for a conservative estimate to the lowest mass halo that retains gas, and so wish to retain spherical symmetry, in order 
to maximises the potential interaction between the host ISM and the SNe. Secondly it is less clear that these very 
small halos include significant angular momentum, and so a non- rotating model is the simplest applicable one 
Taking the mean global temperature of the ISM to be T and the line-of-sight turbulent velocity dispersion to be at, 
the mean density of gas in the potential, <F, is given by (Sutherland & Bicknell 2007): 

/3 

(B5) 

where we take e = 0, and a^ = a^ + kT/[im and gd is the velocity dispersion of the dark matter, assuming also that 

T is the mean ISM temperature, and at is constant throughout the ISM. 

This mean global density distribution is multiplied by the fractal generated above, giving a clumpy distribution with 
the same mean global properties, and which is stable enough over the timescale of the simulations to remain static in 
the absence of a supernova explosion. 


PV) 

p(0) 


= exp 
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Fig. 21.— M60 Off-centred Clumpy explosion commencing at the end of the progenitor ionization phase; see Fig. 20 for details. The star 
is placed at the half gas-mass radius to reflect the fact that most stars form off-centre. The conic ‘zone of avoidance’ sets up in a region 
that is diametrically opposite the progenitor star; a low density nebula extends across most of the core region before the supernova ignites. 
The SN explodes into an off-centred, a warm low-density nebula, superheats the gas through shock heating, and removes much of the gas 
asymmetrically from the dark matter halo out to ryir after 30 Myr. The early zone of avoidance leads to far-side gas re-accreting into the 
core region after the SN event has passed, but less than 10% of the original gas remains within rvir on long timescales. 
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Fig. 22. — M65 Central Clumpy explosion commencing at the end of the progenitor ionization phase. For the uppermost 4 panels (LOO), 
the inner circle has radius rg; the outer circle has radius rvir- Each panel has been magnified in the next 4 panels (LOl); the circle has 
radius rg. Unlike the M60 case, little radiation leaks out from the centre for the first 6 Myr. The ionization is almost entirely restricted to 
the core region (r <rg). The SN explodes into a warm, moderately dense nebula, superheats the gas through shock heating, but the gas 
cools as the nebula expands, which checks its progress. A diffuse ionised nebula forms that extends across the halo but is still bound by it. 
The core region cools and condenses, leading to re-accretion after the bubble collapses starting as early as 3 Myr after the SN explosion. 
All images are central plane density slices at z = 0.0. Log of the normalised density variable in units of 1.0 x 10“^^ g cm“^: greyscale 
minimum = —3.0 (white), maximum = 2.0 (black). A scale bar in parsecs is shown. 
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Fig. 23. — M65 Off-centred Clumpy explosion commencing at the end of the progenitor ionization phase; see Fig. 22 for details. The 
star is placed at the half gas-mass radius to reflect the fact that most stars form off-centre. The SN explodes into an off-centred, warm 
moderately dense nebula and superheats the gas through shock heating. This sends a blast wave into the surrounding nebula but fails to 
maintain the hot cavity that starts to collapse after 6 Myr. The early zone of avoidance leads to far-side gas re-accreting into the core 
region after the SN event has passed. On 30 Myr timescales, there is essentially no gas loss from the dark halo. 
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Fig. 24.— M70 Central Clumpy explosion commencing at the end of the progenitor ionization phase. For the uppermost 4 panels (LOO), 
the inner circle has radius Ts] the outer circle has radius ryir- Each panel has been magnified in the next 4 panels (LOl); the circle has 
radius rg. Much like the M65 case, little radiation leaks out from the centre for the first 6 Myr. The ionization is almost entirely restricted 
to the core region (r <rs). The SN explodes into a warm, moderately dense nebula, superheats the gas through shock heating, but the gas 
cools as the nebula expands, which checks its progress. A diffuse ionised nebula forms that extends across the halo but is still bound by 
it. The core regions cool and condense, leading to re-accretion after the bubble collapses starting as early as 3 Myr after the SN explosion. 
All images are central plane density slices at 2 : = 0.0. Log of the normalised density variable in units of 1.0 X 10“^^ g cm“^: greyscale 
minimum = —3.0 (white), maximum = 2.0 (black). A scale bar in parsecs is shown. 
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Fig. 25. — M70 Off-centred Clumpy explosion commencing at the end of the progenitor ionization phase; see Fig. 24 for details. After 
the progenitor phase, only about half the gas (i.e. one hemisphere) is ionised due to shadowing from the core region which becomes denser 
due to compression by the overpressured ionised volume. The star is placed at the half gas-mass radius to reflect the fact that most stars 
form off-centre. The SN explodes into an off-centred, warm moderately dense nebula and superheats the gas through shock heating. This 
sends a blast wave into the surrounding nebula but fails to maintain the hot cavity that starts to collapse after 6 Myr. The early zone of 
avoidance leads to far-side gas re-accreting into the core region after the SN event has passed. On 30 Myr timescales, there is essentially 
no gas loss from the dark halo. 














